Optimized color conversion

ABSTRACT

Methods and systems herein provide for optimized color conversion. Such color conversion includes identifying color regions in the image data and identifying portions of the color regions that overlap. The color conversion also includes converting color values in the color regions to the color space of the output imaging device based on a plurality of numerical models. The color conversion further includes processing the overlapping portions of the color regions to remove abrupt transitions between the overlapping portions of the color regions. The numerical models may be generated based on human perception and a source of printer errors.

BACKGROUND OF THE INVENTION

1. Field of the Invention

The invention generally relates to the field of printing systems, and inparticular, to methods and systems for performing color conversion.

2. Statement of the Problem

In color printing, displaying, and reproduction, the term gamutrepresents the set of colors that a color-reproduction device isphysically able to generate. Every device that displays or reproduces animage, such as a printer, monitor, scanner, or digital camera, may haveits own unique color gamut. When an image is transferred from one deviceto another, the color gamut of each device is examined to closely matchthe color in the transferred image. That is, an attempt is made toclosely match the color gamut of the image originating device in thedevice to which the image is being transferred so as to provide the mostaesthetically pleasing color conversion. For example, the color gamut ofa digital camera is generally greater than the color gamut of a printer.When color values of the digital camera color gamut are mapped to thecolor gamut of the printer, the conversion process generally requiresintense analysis to ensure that the print quality is sufficiently high.

Adding to the complexity are the color spaces themselves. For example,perceptual color spaces, such as CIELab, are visualized as threedimensional color spaces, where every color that humans can see isuniquely located. Though the CIELab color space is a perceptual colorspace, it is not a perceptually uniform color space as the Euclideandistance in the space does not correspond to the perceptual distance.For example, the magnitude of the perceptual color difference generallydepends upon the color location and the changing direction in chroma andhue. CIELab increasingly overstates the magnitudes of perceived chromadifferences. The human visual system (HVS) is sensitive to the change ofchroma in the neutral color area and insensitive to the change of thechroma in a highly saturated color area. The CIELab color space is alsonon-uniform regarding hue angle in that the thresholds of visualtolerances are a function of hue angle. If the non-uniformity of theCIELab color space in chroma and hue is examined from anotherperspective, CIELab colors have different characteristics at differentlocations. The change of these characteristics is generally continuous.

Most printers (i.e., toner and ink based printers) are able to produceonly a limited numbers of gray levels. In order to produce continuoustone imagery that contains an infinite range of colors or grays, areprographic technique called halftoning is applied to create theillusion of continuous tone images through the use of dot arrangementsand dots of varying size. The combination of the printer halftone designand the specific toner/ink selection determine the number of colors thata printer is physically able to produce (i.e., the gamut).

Printer color conversions are generally performed between adevice-dependent color space and a device-independent color space. Forexample, CMYK color patches whose values span the CMYK color space arefirst printed. The printed patches are then measured using a spectralphotometer which determines the spectral reflectance of each patch undera standard illuminant. Software in the spectral photometer calculatesthe tristimulus values and converts these values to CIELab values. TheCIELab values are floating point numbers. Color conversion lookup tablesare generated using this set of data to convert color values of theCIELab space to the CMYK space. However, the relationship between CMYKcolor space and the CIELab color space is highly nonlinear due to theinteractions of cyan, magenta, yellow, and black planes. The colorconversions, therefore, generally require complex functions. Anothercomplication to the printer color conversion regards the halftonedesign. Most halftone techniques are capable of creating color graylevels continuously but have an abrupt change between two adjacentlevels. Thus, the color conversion between CMYK and CIELab may becontinuous but not differentiable.

Moreover, certain CMYK values may be mapped to the same CIELab value.Since the range of CIELab color space is generally much larger than therange of the CMYK color space, there are no CMYK conversions for manyCIELab values. Additionally, the CIELab to CMYK conversion has onedegree of freedom such that one CIELab value may be mapped to severalCMYK values. Since no function exists for the CIELab to CMYK colorconversion, the conversion is generally achieved by converting theCIELab values to CMY values first with K, or black toner, being lateradded to each CMY value, yielding inaccurate color conversions that arenot aesthetically pleasing. In any case, prior color conversiontechniques do not consider errors associated with printer sources orhuman perceptual tolerances, each of which plays an important role incolor conversion accuracy. Accordingly, there exists a need to convertcolors while ensuring that the converted colors are aestheticallypleasing and vivid.

SUMMARY OF THE SOLUTION

Embodiments of the invention solve the above and other related problemswith methods and associated systems by providing optimized colorconversions. In one embodiment, a method of converting image data to acolor space of an output imaging device includes identifying colorregions in the image data and identifying portions of the color regionsthat overlap in the image data. The method also includes processing theoverlapping portions of the color regions to remove abrupt transitionsbetween the overlapping portions of the color regions and convertingcolor values in the color regions to the color space of the outputimaging device based on a plurality of numerical models. The numericalmodels may be generated according to a plurality of radial basisfunctions for error approximation of each color region in the imagedata, including human perception and a source of printer errors. Forexample, each color region may include L₁ error approximations forneutral color regions, L₂ error approximations for medium range colorregions, and L_(∞) error approximations for saturated color regions.Processing the overlapping portions may entail optimizing theoverlapping portions by converting color values of one color region tothe color space using color values of an adjacent region. Further, themethod may include validating color conversion of the numerical modelsprior to converting the color values to select color conversionalgorithms for the color regions.

In another embodiment, a color conversion system includes a modelingmodule operable to identify color regions in image data and determineportions of the color regions that overlap. The modeling module isfurther operable to generate a plurality of numerical models used inconverting the color regions to a color space of an output imagingdevice. The system also includes a converter operable to convert colorvalues in the color regions to the color space of the output imagingdevice based on the numerical models, wherein color conversion of thecolor values of one color region is performed according to a colorconversion algorithm for the one color region that includes color valuesfrom an overlapping color region.

In another embodiment, a printing system includes a printer and aprinter controller operable to process image data having a plurality ofcolor values and convert the color values of the image data to a colorspace of the printer. The printer controller is further operable toidentify color regions in the image data, determine portions of thecolor regions that overlap, generate a plurality of numerical modelsused in converting the color regions to the color space of the printer,and convert the color values of the color regions to the color spacebased on the numerical models for printing to the printer. The colorconversion of the color values of one color region is performedaccording to a color conversion algorithm for the one color region thatincludes color values from an overlapping color region.

The invention may include other exemplary embodiments described below.

DESCRIPTION OF THE DRAWINGS

The same reference number represents the same element or the same typeof element on all drawings.

FIG. 1 is a block diagram of a color conversion system in one exemplaryembodiment of the invention.

FIG. 2 illustrates a printing system employing a color conversion systemin one exemplary embodiment of the invention.

FIG. 3 illustrates the non uniformity of CIELab chroma.

FIG. 4 is a graphical illustration of the Minimax CharacterizationTheorem (MCT).

FIG. 5 is a graph of exemplary toner/ink values used in color conversionin one exemplary embodiment of the invention.

FIG. 6 is a flowchart illustrating a method of optimizing colorconversion in one exemplary embodiment of the invention.

DETAILED DESCRIPTION OF THE INVENTION

FIGS. 1-6 and the following description depict specific exemplaryembodiments of the present invention to teach those skilled in the arthow to make and use the invention. For the purpose of this teaching,some conventional aspects of the invention have been simplified oromitted. Those skilled in the art will appreciate variations from theseembodiments that fall within the scope of the present invention. Thoseskilled in the art will appreciate that the features described below canbe combined in various ways to form multiple variations of the presentinvention. As a result, the invention is not limited to the specificembodiments described below, but only by the claims and theirequivalents.

FIG. 1 is a block diagram of a color conversion system 100 in oneexemplary embodiment of the invention. In this embodiment, the system100 includes a modeling module 101 which is used to identify colorregions in image data and determine portions of the color regions thatoverlap. The modeling module 101 is further operable to generate aplurality of numerical models used in converting the color use to acolor space of an output imaging device, such as a printer, a computermonitor, or the like. The numerical models may be generated using radialbasis functions, the details of which are described below. The colorconversion system 100 also includes a converter 102 that is operable toconvert color values in the color regions to the color space based onthe numerical models. In doing so, the converter 102 converts the colorvalues of one color region according to a color conversion model for thecolor region by including color values from an overlapping region.

The converter 102 may optimize printer color conversion based on thelocation of colors in a CIELab color space as well as printer sources oferror and human perceptual tolerance. For example, the modeling module101 may approximate neutral, medium, and saturated color regions. Theconverter 102 optimizes the boundaries of adjacent such color regions toensure continuity and smoothness of color conversions of and between thecolor regions. The converter 102 then converts the color values of theimage data to the desired color space.

In generating the numerical models of the color regions, the modelingmodule 101 takes into consideration certain error sources that may beinvolved in the color conversion process. Generally, these errors aremeasured via a Euclidean distance in a CIELab color space, dE, andinclude, for example, measurement errors, machine calibration errors,errors due to printer paper uniformity (e.g., nonuniform fibers, paperroughness, inconsistent paper fillers, etc.), and numerical modelerrors.

Measurement errors and printer paper uniformity errors are generallysmaller than 1 dE. Errors in well calibrated printers are often smallerthan 2 dE, while errors in color conversion models can vary greatly. Forexample, inside gamut color conversion typically has averageinterpolation errors of about 3 dE to 4 dE but is generally capable ofincreasing between 10 dE to 20 dE.

Human perceptual tolerance varies from location to location in theCIELab color space. The modeling module 101 takes such intoconsideration for each of the color regions (i.e., the neutral, medium,and saturated color regions). For color values in the neutral colorregions, the visual color difference tolerance is generally low(approximately 1 dE) and the interpolation errors are generally small.If the accumulated errors from measurement errors, printer paperuniformity errors, and machine calibration errors are added to themodeling data, the effect of the outliers in the modeling is minimized.For the medium color values of increasing chroma values, the perceptualtolerance increases. Outliers still exist, but they play a lessimportant role in the color conversion and generally require anapproximation model that reduces the least squares error (e.g., a normalequation, a singular value decomposition, or a QR factorization). Forcolors in the highly saturated color regions, the visual colordifference tolerance is much higher. The effect of measurement errors,printer paper uniformity errors, and machine calibration errors do notsignificantly contribute to color differences. Accordingly, an RBF basedapproximation may also be used to reduce the maximum error in theinterpolation.

In general, color conversion approximation algorithms of the pastattempted to minimize least square errors. However, those approximationalgorithms still provided inaccurate color conversions. The L_(p) errorapproximations described herein provide more accurate and aestheticallypleasing color conversions. For example, an L₁ error approximation modelprovides for a robust neutral area color conversion whereas an L₂ errorapproximation provides for color conversion in the intermediate regionsand the L_(∞) error approximation provides for the color conversion andthe reduction of the maximum error in the saturated color regions. Whengenerating these numerical models, the modeling module 101 may takeprinter sources of error and human perceptual tolerance intoconsideration.

Continuity and smoothness are generally important to the quality of acolor conversion. The converter 102 assists in this manner by includingdata in adjacent color regions, or clusters, when processing aparticular cluster. In this regard, the converter 102 may optimize theseregions via an optimization algorithm, such as Newton's method, toensure that the overlapping regions are smooth and continuous (i.e.,ensure that the color transitions between color regions do notunnecessarily include abrupt transitions). Exemplary details of thisoperation are described below in the section entitled “FITTINGNON-LINEAR SCATTERED DATA”.

In assisting with the color conversions, the system 100 also includes avalidation module 103. For example, the validation module 103 mayprocess color values of a color region using the numerical modelgenerated by the modeling module 101. In doing so, the validation module103 may generate color conversion estimates based on the color values inthe numerical model. From there, the validation module 103 may optimizecolor conversion estimates generated as L₁ error approximations forneutral color regions, L₂ error approximations for medium range colorregions, and L_(∞) error approximations for saturated color regions, thedetails of which are described below. Based on these optimizations, thevalidation module 103 may determine whether the color conversion isaccurate enough to proceed.

The validation module may also improve accuracy of the numerical modelvia certain validation techniques, such as a cross validation. Forexample, in a “K-fold” cross validation, the validation module 103 mayuse the L₁, L₂, and L_(∞) error approximations to select a K value forthe color conversion. In a CMYK to CIELab conversion, the validationmodule 103 may divide the modeling set into K-fold partition data sets.K-1 folds are then used for training with the remaining one for testingand the K-fold cross validation is applied to determine the best K valueand model. Details of this and other validation techniques are explainedin greater detail below.

FIG. 2 illustrates a printing system 200 employing a color conversionsystem 100 that is used to convert color values in image data to thecolor gamut of a printer, in one exemplary embodiment of the invention.In this regard, the color conversion system 100 may be implemented assoftware instructions that direct the printer controller 202 to performthe color conversion in a manner as described in the system 100 ofFIG. 1. The printing system 200 receives print jobs from one or morehost systems 201 that are communicatively coupled to the printing system200. Thus, when a host system user wishes to print on the two-sidedmedium 213, the host system user may generate a print job that directsthe rasterizer 208 to do so. The printing system 200 includes a printercontroller 202 that receives the print job and rasterizes the print jobvia the rasterizer 208. Rasterization generally regards the conversionof an image described in a vector graphics format into a raster image ofpixels such as a bitmap for output to a video display or printer.Examples of the two-sided medium 213 include paper, transparent film,and any other type of material capable of being imprinted.

Prior to rasterization, the system 100 receives the image data of theprint job and performs a color conversion of the color values thereinsuch that they are aesthetically presentable on the two-sided medium213. As mentioned, the system 100 may direct the printer controller 202to process image data having a plurality of color values and generate anumerical model to convert the color values of the image data to a colorspace of the printer. The printer controller 202 then identifies colorregions in the image data and determines portions of the color regionsthat overlap. Afterwards, the printer controller 202 generates aplurality of numerical models used in converting the color regions tothe color space and then converts the color values in the color regionsto the color space based on the numerical models for printing to theprinter. The color conversion of the color values of one color region isgenerally performed according to a color conversion algorithm for onecolor region by including color values from an overlapping color region.These converted color values may then be optimized to ensure smoothnessand continuity between the color regions. For example, the converter 102may remove abrupt color transitions between overlapping color regionsprior to color conversion. Those skilled in the art should readilyrecognize that the printer controller 201 may be implemented insoftware, hardware, firmware, or a combination thereof to perform thevarious functions associated with the printer controller 202 includingthe novel aspects of color conversion as described herein.

Before discussing the numerical models for CIELab to CMYK and CMYK toCIELab color conversions in detail, a discussion of color spaces is nowpresented to guide the reader through an understanding of colorconversions and gamut mapping. Perceptual color spaces, such as CIEXYZand CIELab, are visualized as three dimensional color spaces, whereevery color that humans can see is uniquely located. The location of anycolor in the space is determined by its color coordinates. In the CIEXYZcolor space, the components are X, Y, and Z, where Y is luminance, X andZ do not correlate to color appearances. In the CIELab color space: thecomponents are L*, a, and b*, where L* is lightness of the color (e, g.,L*=0 yields black and L*=100 indicates white), a* is the value betweenred and green (e.g., negative values indicate green and positive valuesindicate red), and b* is the value between yellow and blue (e.g.,negative values indicate blue and positive values indicate yellow). TheCIELab color space can also be described using cylindrical coordinatesin the L*C*h color space where C*_(ab) is the chroma coordinate, and his the hue coordinate. Below are the mathematical expressions to converta* and b* into C*_(ab) and h:

$\begin{matrix}{C_{ab}^{*} = \sqrt{a^{*2} + b^{*2}}} & \left( {{Eq}.\mspace{14mu} 1} \right) \\{h = {\tan \left( \frac{b^{*}}{a^{*}} \right)}} & \left( {{Eq}.\mspace{14mu} 2} \right)\end{matrix}$

The relationship between CIEXYZ and CIELab may be described using thefollowing equations:

$\begin{matrix}{{L^{*} = {{116\; {f\left( \frac{Y}{Y_{n}} \right)}} - 16}},} & \left( {{Eq}.\mspace{14mu} 3} \right) \\{{a^{*} = {500\left( {{f\left( \frac{X}{X_{n}} \right)} - {f\left( \frac{Y}{Y_{n}} \right)}} \right)}},} & \left( {{Eq}.\mspace{14mu} 4} \right) \\{{b^{*} = {200\left( {{f\left( \frac{Y}{Y_{n}} \right)} - {f\left( \frac{Z}{Z_{n}} \right)}} \right)}},} & \left( {{Eq}.\mspace{14mu} 5} \right)\end{matrix}$

where X_(n), Y_(n), Z_(n) are the tristimulus values of a referencewhite and generally referred to as the white point, described as

$\begin{matrix}{{f(x)} = \left\{ \begin{matrix}x^{\frac{1}{3}} & {x > 0.008856} \\{{7.787x} + \frac{16}{116}} & {x \leq {0.008856.}}\end{matrix} \right.} & \left( {{Eq}.\mspace{14mu} 6} \right)\end{matrix}$

The color difference in either CIEXYZ or CIELab is generally calculatedusing the Euclidean distance. Though the CIEXYZ and CIELab color spaceare perceptual color spaces, they are not perceptually uniform colorspaces as the Euclidean distance in the space does not correspond toperceptual distance. For example, the colors that have the sameEuclidean distance with respect to an anchor color are on a circle.However, these circles are transformed to ellipsoids when the humanperceptual distances are measured. The major axis of the ellipsoidrepresents the higher tolerance of the color difference and the minoraxis of the ellipsoid represents the lower tolerance for the colordifference. The CIELab color space is a perceptually more uniform colorspace than the CIEXYZ color space. Thus, the CIELab color space is thecolor space used most often for color-difference evaluations.

In the CIELab color space, the magnitude of the perceptual colordifference generally depends upon the color location and the changingdirection in chroma and hue. CIELab increasingly overstates themagnitudes of perceived chroma differences. For example, the humanvisual system (HVS) is sensitive to the change of chroma in the neutralcolor area. The HVS is insensitive to the change of the chroma in ahighly saturated color area. FIG. 3 illustrates such change with theellipsoid size increases 301.

The CIELab color space is also non-uniform regarding hue angle in thatthe thresholds of visual tolerances are a function of hue angle. If thenon-uniformity of the CIELab color space in chroma and hue is examinedfrom another perspective, CIELab colors have different characteristicsat different locations. The change of these characteristics is generallycontinuous. Color difference models correlate the perceptual distancewith the Euclidean distance in the CIElab color space based on vastamounts of color difference experimental data.

There are several CIE perceptual-based color difference models. Forexample, CIEDE2000 calculates the perceptual difference by weightedlightness difference, chroma difference, hue difference, and theinteraction between the chroma and the hue difference. However, theCIEDE2000 was developed for small color differences and does not performwell for large color differences. The empirical weighted colordifferences of lightness, chroma, and hue are often used for CIEDE2000,with higher weighting factors for hue and lightness differences andlower weighting factor for chroma differences. The color differencecalculated using the Euclidean distance in CIELab is called deltaE,(also referred as ΔE, or dE_(ab).) dE₀₀ represents the color differencecalculated using CIEDE2000.

In general, there are two types of color management systems among colordevices: additive systems such as a computer monitor and subtractivesystems such as a printer. In additive systems, colors are transmittedlights. Red, green and blue lights are referred as the additive primarycolors. When used in various degrees of intensity and variation, theycreate all other colors of light. When superimposed equally, they creategray. Primary colors of subtractive systems like printers may be cyan,magenta, yellow, and black that are used together to effectively createa multitude of other colors based on the subtractive color theory.

The color conversions for an additive system are normally between theRGB color space and CIEXYZ color space. These conversions are relativelystraightforward since the gamma corrected RGB primaries are linearlyrelated to XYZ values in the CIEXYZ space for most devices. An exampleof such conversion is the conversion for the RGB colors described asfollows:

$\begin{matrix}{{\begin{pmatrix}X \\Y \\Z\end{pmatrix} = {\begin{pmatrix}0.4124 & 0.3576 & 0.1805 \\0.2126 & 0.7152 & 0.0722 \\0.0193 & 0.1192 & 0.9505\end{pmatrix}\begin{pmatrix}{g\left( R_{srgb} \right)} \\{g\left( G_{srgb} \right)} \\{g\left( B_{srgb} \right)}\end{pmatrix}}},{where}} & \left( {{Eq}.\mspace{14mu} 7} \right) \\{{g(k)} = \left\{ {{{\begin{matrix}{\left( \frac{K + \alpha}{1 + \alpha} \right)^{\gamma},} & {K > 0.4045} \\{\frac{K}{12.92},} & {{otherwise},}\end{matrix}\gamma} = 2.2},{{{and}\mspace{14mu} \alpha} = {0.055.}}} \right.} & \left( {{Eq}.\mspace{14mu} 8} \right)\end{matrix}$

As is known, color printers may be toner-based or ink-based. Most of theprinters are only able to produce a limited numbers of gray levels(e.g., a binary printer produces two levels (i.e., 0 and 1) 2-bitprinters produce levels of 0 1, 2, and 3). In order to producecontinuous tone imagery that contains an infinite range of colors orgrays, a reprographic technique called halftoning is applied to createthe illusion of continuous tone images through the use of dotarrangements and dots of varying size. The combination of the printerhalftone design and the specific toner/ink selection determine thenumber of colors that a printer is physically able to produce which iscalled its gamut. The gamut is measured in the CIELab color space.

Printer color conversions depend the toner/ink selections and thehalftone design. The conversions are generally performed between adevice-dependent color space (CMYK) and a device-independent color space(CIELab). Often the format of the color conversion is an equally spacedlookup table, the values of which are either 1-byte or 2-byte integers.

In generating a printer color conversion, CMYK color patches whosevalues span the CMYK color spaces are printed. Each component value ofCMYK is normally a 1-byte integer. The printed patches are measuredusing a spectral photometer which measures the spectral reflectance ofeach patch under a standard illuminant. The software in the spectralphotometer calculates the tristimulus values and converts these valuesto CIELab values. The CIELab values are floating point numbers. Followedby the implementation of a data modeling algorithm and a gamut mappingalgorithm, the color conversions LUTs are generated using this set ofdata. The CMYK and the CIELab values are normalized to the range of 0 to1 for these calculations. Both CMYK to CIELab and CIELab to CMYKconversions may be generated for a color printer. The CMYK to CIELabcolor conversion is used when an image generated in this CMYK colorspace is to display on a monitor or to print on a different colorprinter. In other words, this conversion is used when the printer CMYKcolor space is used as a source color space. The CIELab to CMYK colorconversion is used when the CIELab color space is the source colorspace.

For many color printers, the relationship between CMYK color space andthe CIELab color space is highly nonlinear due to the interactions ofcyan, magenta, yellow, and black planes in the printer subtractivesystem, thereby requiring complex functions. Another complication to theprinter color conversion is due to the halftone design. Most halftonetechniques are capable of creating color gray levels continuously.However, an abrupt change of two adjacent levels is a weakness in manyhalftone designs. Thus, the color conversion between CMYK and CIELab maybe continuous but not differentiable.

The CMYK to CIELab color conversion is neither injective nor subjective.For example, every color in the CMYK color space may be mapped to acolor in the CIELab color space. However, different CMYK values may alsobe mapped to the same CIELab value. So the CMYK to CIELab colorconversion is not injective. The range of CIELab color space is muchlarger than the range of CMYK color space. That is, for many CIELabvalues, there are no CMYK values such that f(C, M, Y, K)→(L, a, b). Thusthe CMYK to CIELab color conversion is not subjective.

Additionally, the CIELab to CMYK conversion is from R3 to R4. Thus, onlyone degree of freedom in this conversion exists and one CIELab value maybe mapped to several CMYK values. On the other hand, because of thesmaller gamut of the color printer, not every CIELab value can be mappedto a CMYK value. Since no function exists for the CIELab to CMYK colorconversion, the conversion is achieved by converting the CIELab valuesto CMY values first and a special rule is created for adding the blacktoner to each CMY value.

The CIELab to CMY color conversion within a printer gamut is bothinjective and subjective. That is, virtually every CIELab value insidethe printer color gamut can be mapped to a unique color in the CMY colorspace. So the CIELab to CMY color conversion is injective within colorgamut. Additionally, for CMY values, there exists a CIELab value suchthat g (L, a, b)→(C, M, Y) where the CIELab to CMY color conversion issubjective. Out-of gamut CIELab values generally cannot be mapped to aCMY value. Therefore, no functions exist for the CIELab to CMYconversion outside the printer gamut and color conversion between theCMYK and CIELab color spaces is achieved by approximation.

Understanding the sources of error is important for the design of anapproximation model. In addition, it is important to understand howthese errors interfere with human perceptions for the design of theprinter color conversion model. There are four significant sources oferror in the color conversion process. These errors are normallymeasured with the Euclidean distance in the CIELab color space, ΔE, andinclude:

-   -   1. measurement errors from the instrument, which are normally        smaller than 1 ΔE;    -   2. errors due to the machine reliability and repeatability,        which for a well-calibrated printer are generally smaller than 2        units of ΔE;    -   3. errors due to paper roughness and non-uniform paper fibers        and fillers, which are typically smaller than 1 unit of ΔE; and    -   4. errors from the model of the color conversion, which average        around 3-4 units of ΔE with a maximum error between about 10 to        20 units of ΔE.

Because human perceptual tolerance varies from location to location inthe CIELab color space, care is taken in the approximation algorithms indifferent color areas. For colors in neutral areas, the visual colordifference tolerance is low (around 1 ΔE). The interpolation errors inthese regions are generally small. If the accumulated errors from thefirst three categories are added to the modeling data, the effect of theoutliers in the modeling set should be minimized. For the medium colors,with the increase of the chroma values in these areas, the perceptualtolerance increases. Though the outliers still exist, they play a muchless important role in the color conversion. It is appropriate to designan approximation model to minimize the least squares error. For colorsin the highly saturated color regions, interpolation errors arerelatively high. The approximation algorithm is used to minimize themaximum error in the interpolation.

In general, color conversion approximation algorithms aim at minimizingleast square errors. Due to the sources of error described above,emphasis is placed on the L_(p) error criteria, in particular the L₁ andL_(∞) error for a robust neutral area color conversion and the reductionof maximum error in the saturated color regions. One approximation modelincludes the natural neighbor interpolation method, which is based onthe Voronoi region. A Voronoi region V_(i) is a set of points for whicha center c_(i) is the nearest center vector, as follows:

V _(i) ={x ε R ^(n) :i=arg min_(jε1) ∥x−c _(j)∥},   (Eq. 9)

where I is the set of center indices. This interpolation method is alocal region-based weighting method generally provides better underlyingprinter characteristics.

Global interpolation methods may also include the inverse distanceweighting method and the RBF method. For the inverse distance weightingmethod, data points are used to interpolate a value by assuming thatinterpolated values are affected more by nearby points and less by moredistant points. The assumption is not accurate for colors in manyregions. Thus, this method has yielded relatively poor color conversionresults. RBF data fitting equipped with an L_(p) error-criterion, on theother hand, improves color conversion. The embodiments described hereinalso take into account perceptual color difference tolerance. Asmentioned, perceptual color difference tolerances generally vary overdifferent CIELab color regions. Taking this into account, the RBF datafitting technique has the advantage of minimizing color differences indifferent areas with respect to human perceptions or yield results.

Although black toner/ink is one of the primaries in a color printer, theblack toner/ink actually serves as a “secondary” toner in colorprinting. Because the mix of cyan, magenta, and yellow toner/inkproduces a “muddy brown” instead of black, black toner/ink is added into make a more rich black and thereby make other colors richer anddarker. Black toner/ink also typically costs less than cyan, magenta,and yellow toner/ink making the embodiments herein more advantageousfrom a cost perspective.

Printer CIELab to CMYK color conversion maps from R3 to R4.Additionally, CIELab is mapped to the CMY color space using aninterpolation technique, then black toner/ink denoted by K, is added.The amount of cyan, magenta, yellow toner/ink may be adjusted byapplying the under color removal (UCR) and the under color addition(UCA) corrections. The typical way of applying UCR is to remove equalamounts of yellow, magenta, and cyan and replace them with the sameamount of black. Gray component replacement (GCR) is the result of UCRfollowed by UCA. However, colors resulting from UCR are generally dark,less saturated, and hue shifted. UCA may be applied to improve theweakness of UCR. Small portions of the three colors of cyan, magenta,and yellow may be added back to these areas. The UCA process isempirical and relatively processor intensive due to the mapping from alower dimensional space to a higher dimensional space (e.g., R3 to R4)and the lack of a high fidelity model describing the nonlinear behaviorof each toner. UCA is normally achieved by printing many patches andchoosing the best results. In practice, the UCR and UCA correctionschemes introduce errors in the color conversion process.

Because the CIELab to CMYK conversion is a map from R3 to R4, a set ofCMYK values corresponding to the same CIELab value may be established.In this regard, g denotes a mapping from the CIELab values to the CMYKvalues such that

∃{(C, M, Y, K)^((i)) ,i=1, 2, . . . . k},s.t.,g ⁻¹(C, M, Y,K)^((i))=(L,a,b),   (Eq. 10)

where, i=1; 2; . . . ; k. These level sets of CMYK values may be foundthrough an optimization process if approximation models are developedfor the CMYK to CIELab conversion and the CIELab to CMY conversion.

As mentioned, gamut mismatch is problematic in color reproductionsystems since each device has a different gamut. In general inputdevices such as digital cameras, scanners, etc., generate more colorsthan the output devices such as color printers. When the input colorspace is larger than the gamut of the output color device, gamut mappingalgorithms are applied. The gamut mapping process transforms a point inthe source gamut to a realizable color inside the gamut of the outputdevice (e.g., a destination gamut). The form of this transformation canimpact the quality of the reproduced images, especially print images.For high-end digital printers, the gamut is relatively small. Themaximum error of color conversions normally derives from mapping anoutside-gamut color to an inside-gamut color mapping.

In the gamut mapping process, a device color gamut boundary may beinitially identified. The calculated boundaries are normally connectedlines. Previous gamut mapping algorithms mapped out-of-gamut pointsdirectly to the destination gamut. The most typical gamut mappingsinclude clipping and scaling algorithms. Clipping algorithms clipout-of-gamut points to the destination gamut boundary. Scalingalgorithms, on the other hand, scale the input color gamut to outputcolor gamut such that some out-of-gamut points are mapped inside thedestination gamut while others are mapped to the boundary of thedestination gamut. The mapping directions may include a choice of linesalong which the mapping is applied (e.g., a line changing in chroma atconstant lightness and hue) that are often decided based on visualexperiments. The mapping directions may also be along the lines ofconstant perceptual attribute predictors (e.g., constant lightness andhue, constant saturation and hue), along the lines towards a singlecenter-of-gravity (e.g., the center is the compression towards L*=50 onlightness axis), along the lines towards variable center-of-gravity(e.g., the cusp is the compression towards lightness of cusp onlightness axis), or along the lines towards the nearest color in thereproduction gamut (e.g. as is the case with minimum ΔE clipping).

The large variability in past color gamut mapping suggests that nouniversal gamut mapping algorithm works for every device and image.Ideal gamut mapping generally depends on image content, preservation ofperceived hue throughout color space, and the extent of the gamutmismatch in various regions of color space. Relationships among colordata developed by a numerical model have not previously been consideredin the gamut mapping process. The underlying relationship among colordata may more accurately retain or rationally transform after gamutmapping. Such is achieved herein using a numerical model with guidancefrom the perceptual model.

The numerical model quantitatively compares the color difference betweenthe original color and the printed color. The original colors are theCIELab colors. These colors are converted with the printer colorconversion to the CMYK colors. The CMYK colors are sent to the printerto print and then measured by spectral photometers. The numerical methodtests the approximation errors inside color gamut and gamut mappingresults out-side the color gamut. The mean error and the maximum errormay be calculated using either dE_(ab) or dE₀₀ for the color conversioninside the printer gamut. The maximum conversion error is often weightedmore than the mean error in the evaluation because the color differenceis more objectionable than the average color difference.

The perceptual model, on the other hand, correlates human perceptionswith physical stimuli. A set of test images including pictorial images,vector graphics, and text samples are printed from the printer. Thecolors in the images include highlight colors, low-light colors,saturated colors, and out-of-gamut colors. Observers look at theseimages and make perceptual judgments based on the selected psychophysicsmethod, such as the paired comparison method. On such experiments,observers are asked to compare an image pair converted by differentcolor conversion methods. The original photos are normally presented asthe reference. Typical questions are “Which one do you like better?”,“Which one is more similar to the original image?”, etc. The data isanalyzed and converted to perceptual values using Thurstone's law ofcomparative judgment.

To assist the reader in the development of the numerical model, thisdisclosure broken into certain sections. Specifically, the sections areentitled “FITTING NON-LINEAR SCATTERED DATA”, “OVERDETERMINED LINEARSYSTEMS IN THE L_(∞) SENSE”, “COLOR CONVERSION USING THE L₁ NORMAPPROXIMATION”, “ONE DIMENSIONAL COLOR MANIFOLDS IN THE CMYK COLORSPACE”, “COLOR CONVERSION WITH INK LIMITATIONS”, and “GAMUT MAPPING”.The section entitled “FITTING NON-LINEAR SCATTERED DATA” provides thereader with a detailed description of solving printer color conversionproblems using RBFs and data fitting techniques adapted to the L_(∞) andthe L₁ error criteria. The section entitled “OVERDETERMINED LINEARSYSTEMS IN THE L_(∞) SENSE” presents the reader with the Barrodale andPhillips dual algorithm (the “BP” algorithm) and the Bartels, Conn, andLi primal algorithm (the “BC” algorithm) for converting relatively largesparse problems to linear programming problems. The section entitled“COLOR CONVERSION USING THE L₁ NORM APPROXIMATION” illustrates a generalL₁ approximation for overdetermined linear systems, the BR algorithm forsolving the primal problem, and the L₁ approximation as it is used incolor conversion. The section entitled “ONE DIMENSIONAL COLOR MANIFOLDSIN THE CMYK COLOR SPACE” introduces the reader to an optimized GCRalgorithm with relatively high quality toner/ink selection. The sectionentitled “COLOR CONVERSION WITH INK LIMITATIONS” discusses toner/inklimitations implemented as a constraint in the linear programmingproblem formulations. The section entitled “GAMUT MAPPING” introducesthe reader to the concept of numerical gamut mapping based on the RBFsdiscussed above and enhanced by the perceptual gamut mapping.

Fitting Non-Linear Scattered Data

Problems with the interpolation techniques are overcome through the useof RBFs that include data fitting techniques adapted to L_(∞) and L₁error criteria. Additionally, a perceptual color model is used to guidethe implementation of the RBFs. For example, in order to solve nonlinearapproximation problem using RBFs, a clustering algorithm is employedusing the combination of the Linde-Buzo-Gray (LBG) algorithm and thecovariance matrix-based algorithm. RBFs may then be selected andparameters thereof may be optimized. Approximation algorithms aregenerated for the L_(p) norm (emphasis on the L_(∞), L₂, and L₁ norm)for the neutral, intermediate, and saturated color areas. A K-fold crossvalidation then refines the model parameters and an optimizedoverlapping area between two different models is determined.

Theoretically, a color printer is capable of producing at least 4 bytesof colors (1 byte each of cyan, magenta, yellow, and black). A subset ofthese CMYK colors and their associated CIELab values may be used togenerate a color conversion model. Because of the underlying complicatedcharacteristics of a color printer, a large CMYK data set and theassociated CIELab values (often above 2000 data points) are selected togenerate the color conversion model. The data selection depends in parton the halftone design which occasionally yields a nonlinear response.The data selection also depends on the perceptual color differencetolerance as the perceptual color difference tolerances vary overdifferent CIELab color regions. Therefore, different color conversionmodels should be generated for different regions based on the perceptualcolor difference tolerance. The continuity and smoothness of these colorconversion models are important to the overall quality of the colorconversion. RBFs provide the flexibility of applying to both theinterpolation and approximation data fitting problems. RBFs are usefulto interpolate scattered data as the associated system of linearequations is invertible under relatively mild conditions on thelocations of the data points.

Generally, RBFs are useful for approximating an unknown function fromthe data in the form

$\begin{matrix}{{{f(x)} = {{Ax} + a_{0} + {\sum\limits_{i = 1}^{N}{w_{i}{\varphi \left( {{x - c_{i}}} \right)}}}}},} & \left( {{Eq}.\mspace{14mu} 11} \right)\end{matrix}$

where f(x) is assumed to be continuous and differentiable over a compactdomain of x ε R^(n) and f(x) ε R The term Ax represents the linearcomponent of the mapping and the offset a₀ allows for the mean of thedata to be non-zero. The nonlinear portion of the map is represented bythe superposition of the nonlinear function φ weighted by vectors w_(i).

The RBF has several important features. For example, the weightparameters may be determined using linear models, the RBFs may be eitherlocal or global and the locations of the basis functions may be adaptedusing different clustering routines. For color conversion, the Voronoiregion, for example, may provide relatively good clustering because itcan incorporate the color difference model to characterize each clusterwith its unique characteristics. Also, RBFs generally do not requiredata to lie on any sort of regular grid and the function f(x) iscontinuous and differentiable.

Traditionally, the weights w_(i) were determined by minimizing themean-square error. However, due to the unique error characteristics incolor conversion, the algorithms for solving the RBFs' weights usingL_(p) norm emphasizes L_(∞) norm, L₂ norm, and L₁ norm. Accordingly, theapproach described herein includes choosing a training/testing data setfor color conversion, clustering the color data, determining the centersoptimizing the number of centers, and optimizing RBFs and parameterselections. Then, the L_(p) norm selections (specifically, the L_(∞),L₂, and L₁ norms) may be optimized and the RBFs combined with linearprogramming for the L_(∞) norm and L₁ norm approximation employing theduality theorem. Thereafter, the perceptual color model may be employedto improve the qualitative result.

The color conversion model is built based on the availability of a dataset of input-output pairs, such as CMYK vs. CIELab. There are many waysto select a set of training data. For example, the CMYK colors inregular grids may be selected or the CMYK colors based on the halftonecurves may be selected. In any case, the CMYK data is printed on a colorprinter and measured in the CIELab color space. For the CMYK to CIELabconversion, these CIELab values are random within the printer gamut andcan be expressed for f: R4→R3 in the following form:

{x ^((i))=(C, M, Y, K)^((i)) ε R4, i=1,2, . . . , m}  (Eq. 12)

{y ^((i))=(L,a,b)^((i)) ε R3, i=1,2, . . . , m},   (Eq. 13)

where m is the number of color pairs in the data set. For the CIELab toCMYK conversion, g : R4→R3

{x ^((i))=(L,a,b)^((i)) ε R3, i=1,2, . . . , m}  (Eq. 14)

{y ^((i))=(C, M, Y, K) ^((i)) ε R4, i=1,2, . . . , m}.   (Eq. 15)

The printer color conversions described herein include both CMYK toCIELab and the CIELab to CMYK color conversions. The CIELab to CMYKconversion is achieved by converting the CIELab values to the CMY valuesfirst and then adding black to each CMY value via the GCR process.Therefore, the CMYK set should have enough C,M,Y, and K values.Generally, the grid points of cyan, magenta, yellow, and blackcomponents are used to represent the curve change of each componentcreated by the halftone process. This means that the grids point may notbe equally spaced.

Data clustering is used to identify the locations where basis functionsshould be placed in the model. In general, RBF algorithms focus onplacing functions such that the mean square error is minimized. Due tothe uniqueness of perceptual color characteristics, the location of thebasis function should be placed in the local regions which have similarproperties in terms of human visual tolerances. As mentioned, the CMYKvalues are device-specific values and do not directly relate to humanperception. In the construction of the mapping g: R4→R3and f: R4→R3,clustering is performed on the CIELab data so the characteristics of thecolor data can be described by each cluster. The visual tolerances ofcolor data may be determined using a color difference model. The mostrecent color difference model CIEDE2000 is used as the distance metricfor each color relative to the center color because the distances aredesigned to be small.

The clustering algorithm also combines the Voronoi-based clusteringalgorithm with a subset selection algorithm. Because the CIELab data arerandom, a random center selection algorithm may be used as the firststep for clustering/center selections. The clustering algorithms arethen used to identify the clusters and to refine the center of eachcluster. There are several efficient clustering algorithms availablethat may be used including “k-means,” LBG, and topology-preservingmappings. The algorithm generally proceeds by selecting an initialnumber N_(c) of centers randomly from the m data points in the trainingset and, for each center, computing the Voronoi set (i.e., all points inthe data closest in the perceptual color space to this center relativeto other centers). Then, the centers as the mean of the points in theVoronoi set are updated. This is repeated until the CIELab values ofthese centers no longer change. The clustering algorithm described abovecan be effectively combined with a subset selection algorithm, such asthe orthogonal least squares (OLS) algorithm or the covariance matrixalgorithm of commonly owned and co pending U.S. Patent ApplicationPublication No. 2007/009796, the entire contents of which areincorporated by reference, to eliminate unnecessary clusters.

The OLS algorithm identifies which of the N centers are most useful inthe RBF model and indicates which centers may be deleted from the modeland at what expense to accuracy expressed in color difference. Thecomputation of the weight parameters in Eq. 11 is an over-determinedleast squares problem. Thus, a solution is sought to the set ofinconsistent equations,

y=Φw,   (Eq. 16)

where each column of the matrix Φ is associated with a single center andf does not actually reside in the column space of Φ. In general, aquality function is associated with a center. In this case, the value ofa column φ_(i) of Φ is measured to solve the least squares. One measureis the cosine of the angle between f and each column of Φ, such as

$\begin{matrix}{v_{i} = {\left( \frac{f^{T}\varphi^{(i)}}{{f}{\varphi^{(i)}}} \right)^{2}.}} & \left( {{Eq}.\mspace{14mu} 17} \right)\end{matrix}$

In other words, the center is the one for which v_(i) is largest (i.e.,where the angle is smallest). This solution is indexed by i* meaningφ_(i)* is the solution to optimization. Once the column and hence theRBF center of Φ has been determined, the next best center can becomputed by projecting the remaining columns along φ_(i)*. This approachmay be iterated to obtain a reduced subset of centers that hasacceptable modeling accuracy. This algorithm may therefore proceeds asfollows:

$\begin{matrix}{{{{{Let}\mspace{14mu} k} = 1},2,\ldots \mspace{11mu},N_{c}}{{{Step}\mspace{14mu} 1\text{:~~}k} = 1}} & \; \\{i_{1} = {{argmax}\left( \frac{f^{T}\varphi^{(i)}}{{f}{\varphi^{(i)}}} \right)}^{2}} & \left( {{Eq}.\mspace{14mu} 18} \right) \\{q^{(1)} = \varphi^{(i_{1})}} & \left( {{Eq}.\mspace{14mu} 19} \right) \\{{{Step}\mspace{14mu} k} = 1} & \; \\{{\alpha_{jk}^{(i)} = \frac{q_{j}^{T}\varphi_{i}}{q_{j}^{T}q_{j}}},\mspace{14mu} {1 \leq j < k}} & \left( {{Eq}.\mspace{14mu} 20} \right) \\{q_{k}^{(i)} = {\varphi_{i} - {\sum\limits_{j = 1}^{k - 1}{\alpha_{jk}^{(i)}q_{j}}}}} & \left( {{Eq}.\mspace{14mu} 21} \right) \\{{i_{k} = {{argmax}\left( \frac{f^{T}q_{k}^{(i)}}{{f}{q^{(i)}}} \right)}^{2}}{{{{where}\mspace{14mu} i_{k}} \leq i \leq N_{c}},{i \neq i_{1}},{{\ldots \mspace{14mu} i} \neq i_{k - 1}},{and}}} & \left( {{Eq}.\mspace{14mu} 22} \right) \\{q_{k} = {q_{k}^{i_{k}} = {\varphi_{i_{k}} - {\sum\limits_{j = 1}^{k - 1}{\alpha_{jk}{q_{k}.}}}}}} & \left( {{Eq}.\mspace{14mu} 23} \right)\end{matrix}$

The algorithm ends at the Mth step when

$\begin{matrix}{{{1 - {\sum\limits_{j = 1}^{M}\left( \left( \frac{f^{T}q_{j}}{{f}{q_{j}}} \right)^{2} \right)}} < \rho},{{{where}\mspace{14mu} 0} < \rho < 1}} & \left( {{Eq}.\mspace{14mu} 24} \right)\end{matrix}$

A cluster can also be described by its mean vector andvariance-covariance matrix. The mean vector is the center of the clusterand the variance-covariance matrix consists of the variances of thecolor data points along the main diagonal and the covariance betweeneach pair of color components in the other matrix positions. Forexample, if a cluster is in the CMYK color space, the covariance of thiscluster provides a measure of the strength of the correlation betweenany two color components. Let the center CMYK value of a cluster be x=(C, M, Y, K) (i.e., x _(i)= C, x ₂= M, x ₃= Y. x ₄= K) and let m_(c)denote the number of data points in a cluster. The formula for computingeach element COV_(jj), in a 4×4 CMYK covariance matrix for a cluster is:

$\begin{matrix}{{{COV}_{{jj}^{\prime}} = {\sum\limits_{i = j}^{mc}{\left\lbrack {\left( {x_{j}^{i} - {\overset{\_}{x}}_{j}} \right)\left( {x_{j^{\prime}}^{i} - {\overset{\_}{x}}_{j^{\prime}}} \right)} \right\rbrack/\left( {m_{c} - 1} \right)}}},} & \left( {{Eq}.\mspace{14mu} 25} \right)\end{matrix}$

where j;j′=1, . . . , 4. When the rank of the covariance matrix is lessthan the number of the input color components, the correlations are notall independent of each other. Also, the number of the data points inthe cluster is less than the number of color components.

A Gaussian function which can be expressed as the function of thecovariance matrix is not selected if a non-full rank cluster exists.Thus, the clustering algorithm proceeds by calculating the covariancematrix of each cluster and, if the rank of the covariance matrix of acluster is less than the number of input color components, the center ofthis cluster is deleted and the data points in the cluster are moved tonearby clusters. Thereafter, the center of each cluster is updated andthe LBG algorithm is updated to refine each cluster. This is repeateduntil every cluster has a full rank covariance matrix.

The radial basis functions themselves are generally selected so thatthey satisfy an invertibility condition on the square interpolationproblem for the weights. These functions have expansions, such as thosegiven by Eq. 11, that represent continuous functions over compactdomains. The location of these functions is stipulated by the vectorcenters {ck} generally scattered over the domain in a manner thatreflects the distribution of the data. Both theoretical investigationand practical results suggest that the choice of the nonlinearity φ(.)is not critical to the performance of the RBF. However, for the colorconversions from CMYK to CIELab, and CIELab to CMYK, the multiquadricfunction

φ(x)=√{square root over (r ² +x ²)}  (Eq. 26)

is generally more accurate. Other functions that may be used include theGaussian

$\begin{matrix}{{\varphi (x)} = {\exp \left( \frac{- x^{2}}{r^{2}} \right)}} & \left( {{Eq}.\mspace{14mu} 27} \right)\end{matrix}$

and the inverse multiquadratic function

$\begin{matrix}{{\varphi (x)} = {\frac{1}{\sqrt{r^{2} + x^{2}}}.}} & \left( {{Eq}.\mspace{14mu} 28} \right)\end{matrix}$

Regardless of the function used, the radius r in Gaussian,multiquadratic, and inverse multiquadratic is optimized. Generally, themultidimensional Newton's method is an efficient optimization method for2-norm least square settings. In this case, the objective function aimsat minimizing the mean error. For example, using the color conversionfrom CMYK to CIELab:R4→R3,

$\begin{matrix}{{f_{L}\left( {x_{i},r_{L}} \right)} = {{A_{L}x_{i}} + a_{L_{0}} + {\sum\limits_{j = 1}^{n}{w_{L_{j}}{\varphi_{L}\left( {{{x_{i} - c_{j}}},r_{L}} \right)}}}}} & \left( {{Eq}.\mspace{14mu} 29} \right) \\{{f_{a}\left( {x_{i},r_{a}} \right)} = {{A_{a}x_{i}} + a_{a\; 0} + {\sum\limits_{j = 1}^{n}{w_{aj}{\varphi_{a}\left( {{{x_{i} - c_{j}}},r_{a}} \right)}}}}} & \left( {{Eq}.\mspace{14mu} 30} \right) \\{{f_{b}\left( {x_{i},r_{b}} \right)} = {{A_{L}x_{i}} + a_{b\; 0} + {\sum\limits_{j = 1}^{n}{w_{bj}{\varphi_{b}\left( {{{x_{i} - c_{j}}},r_{b}} \right)}}}}} & \left( {{Eq}.\mspace{14mu} 31} \right)\end{matrix}$

The problem formulation is thus:

$\begin{matrix}{{E\left( {r_{L},r_{a},r_{b}} \right)} = {\sum\limits_{i = 1}^{m}\left\lbrack {{f_{L}\left( {x_{i},r_{L}} \right)} - {y_{1}^{i}\text{)}^{2}} + {f_{a}\left( {x_{i},r_{a}} \right)} - {y_{2}^{i}\text{)}^{2}} + {f_{b}\left( {x_{i},r_{b}} \right)} - {y_{3}^{i}\text{)}^{2}}} \right\rbrack}} & \left( {{Eq}.\mspace{14mu} 32} \right) \\{E_{\min} = {E\left( {r_{L},r_{a},r_{b}} \right)}} & \left( {{Eq}.\mspace{14mu} 33} \right)\end{matrix}$

Then for each iteration, let

$\begin{matrix}{g^{(k)} = \begin{pmatrix}{{\partial E}/r_{L}} \\{{\partial E}/r_{a}} \\{{\partial E}/r_{b}}\end{pmatrix}} & \left( {{Eq}.\mspace{14mu} 34} \right) \\{{F(k)} = {\begin{pmatrix}\frac{\partial^{2}E}{\partial^{2}r_{L}} & \frac{\partial^{2}E}{{\partial r_{L}}{\partial r_{a}}} & \frac{\partial^{2}E}{{\partial r_{L}}{\partial r_{b}}} \\\frac{\partial^{2}E}{{\partial r_{a}}{\partial r_{L}}} & \frac{\partial^{2}E}{\partial^{2}r_{a}} & \frac{\partial^{2}E}{{\partial r_{a}}{\partial r_{b}}} \\\frac{\partial^{2}E}{{\partial r_{b}}{\partial r_{L}}} & \frac{\partial^{2}E}{{\partial r_{b}}{\partial r_{a}}} & \frac{\partial^{2}E}{\partial^{2}r_{b}}\end{pmatrix}.}} & \left( {{Eq}.\mspace{14mu} 35} \right)\end{matrix}$

Accordingly,

r _(k+1) =r ^((k)) −F(r ^((k)))⁽⁻¹⁾ g ^((k))   (Eq. 36)

For every r^((k+1)), the weighting factor w_(j)'s are updated and thenew function E(r) is obtained. The iteration stops when E(r) does notimprove.

The algorithms for RBF data fitting may depend on the selection ofnorms. For color conversions, the norm may be selected considering thevisual tolerance of the colors in the CIELab color space and the sourceof errors in the color conversion process. For example, different normsshould be selected for differential color regions, such as

L₁ norm:

$\begin{matrix}{{x}_{1} = {\sum\limits_{i = 1}^{m}{x_{i}}}} & \left( {{Eq}.\mspace{14mu} 37} \right)\end{matrix}$

L₂ norm:

$\begin{matrix}{{x}_{2} = \sqrt{\sum\limits_{i = 1}^{m}x_{i}^{2}}} & \left( {{Eq}.\mspace{14mu} 38} \right)\end{matrix}$

L_(∞) norm:

$\begin{matrix}{{x}_{\infty} = {\max\limits_{1 \leq x \leq m}{{x_{i}}.}}} & \left( {{Eq}.\mspace{14mu} 39} \right)\end{matrix}$

The relationship of errors via L₁norm, L₂ norm, and L_(∞) may be statedas:For all e in C[a; b] the inequalities,

∥e∥ ₁≦(b−a)^(1/2) ∥e∥ ₂≦(b−a)∥e∥ _(∞)  (Eq. 40)

hold.

This theorem indicates that it small errors are found via L_(∞) norm,then the errors via L₁ norm and L₂ norm are also small. Thus, [a, b]=[0,1], and the CIELab color space may be heuristically chosen in threeregions:

-   1. Neutral region: 0≦C*_(ab)≦7;-   2. Medium color region: 7<C*_(ab)≦30; and-   3. Saturated color region: C*_(ab)>30,    where C*_(ab) is the chroma of colors in CIELab color space.

For the neutral color region, the visual color difference tolerance isrelatively small. Error from the measurements and machine can affect theaccuracy of color conversion. The L₁ norm is chosen since it is leastsensitive to the outliers. For the medium color range, the colordifference tolerance is larger and the errors from the measurements andmachine play less of a role affecting the quality of the colorconversions. However, control of the color conversion error should beless than a magnitude of 2 so the color difference is less visible. TheL₂ norm is chosen to minimize the average error. For the saturated colorregions, the errors are generally small compared to the color conversionerror of gamut mapping and it is unlikely for the HVS to detect anyerrors with magnitude less than 2. To minimize the maximum error of thecolor conversion, the La norm is generally chosen.

As mentioned, the continuity and smoothness are important to the qualityof the color conversion. When adjacent clusters or two adjacent regionshave different mapping functions, data in the adjacent clusters may beincluded in the current cluster. Overlapping ranges specified as δ's foreach region are optimized. The δ values are the distances calculated inthe CIELab color space. The data regions partitioned for the L₁, L₂, andL_(∞) approximation are:

-   1. L₁ approximation: 0≦C*_(ab)<7+δ_(n). Intended range: 0≦C*_(ab)<7;-   2. L₂ approximation: 7−δ_(ml)<C*_(ab)<30+δ_(mu). Intended range:    7≦C*_(ab)≦30;-   3. L_(∞) approximation: C*_(ab)>30−δ_(s). Intended range:    C*_(ab)≦30,    where δ_(n) is the extended region for the L₁ norm, δ_(ml) and    δ_(mu) are the upper and lower extended regions for the L₂ norm, and    δ_(s) is the extended region for the L_(∞) norm. For the CMYK to    CIELab conversion, the domain space is the undivided CMYK color    space. For the Cielab to CMY conversion, the domain space is the    CIELab color space divided into three regions based on the selected    set of training data described above. Although shown and described    with a particular range, the invention is not intended to be so    limited. For example, the approximation may be varied as desired to    change the granularity of the approximation.

Validation techniques, such as “cross validation” and “bootstraps”, maybe used for model accuracy estimation and model selection. Among thecross validation techniques, the “K-fold cross” validation and the“leave-one-out” cross validation may be used. Using the cross validationtechnique, the absolute accuracy and the bias of variances should bebalanced.

In general, several thousand data points are selected for building acolor conversion model. A relatively large K number may be chosen toensure the model accuracy. Different from other cross validationapplications, the evaluation of color conversion models involves bothnumerical analysis and psychophysics analysis. For the numericalanalysis, the test data set is collected spanning the input color space.For the perceptual analysis, a number of images may be selected in avisual experiment.

L₁, L₂, and L_(∞) approximations may be used to select the K value tothe color conversion. For the L₂ approximation, the mean error may beused as the measure for the model selection. And, for the L₁ and L_(∞)approximations, the L₁ error and maximum error may be calculatedinstead, respectively. In any case, the variance as a function of the Kvalue may be obtained to decide the K value. The cross validationtechnique is described using CMYK to CIELab conversion as an example.First, the modeling set is divided into K-fold partition data sets. Foreach of K experiments, K-1 folds are used for training with theremaining one for testing. Then, K-fold cross validation is applied todetermine the best K value and model. The mean error decreases with theincrease of the K value for K<10 and it stays roughly at the same levelfor K>15. On the other hand, the variance of the mean error staysrelatively flat for K<15 and increases with the increase of the K value.Although the maximum error decreases and the mean error stays relativelyflat for K≦15, the K value for the least squares approximation for thecolor conversion is generally 15 due to the increase of the variance.

The CMYK to CIELab conversion may be performed using the RBF datafitting algorithm in the L₂ norm. In this regard, data may be clusteredusing the combination of the LBG algorithm and the covariance matrixalgorithm. Generally, there are two optimization stages. Initially, arelatively good starting function and radius for the R4→R3 optimizationis sought. Then, models of R4→R for the CMYK to L*, a*, b* conversionsare sought.

In the next modeling stage, f_(L), f_(a), and f_(b) are set tomultiquadratic with initial r_(L), r_(a), and r_(b). Themultidimensional Newton's method is then applied to optimize theseparameters. For the CIELab to CMY conversion, both data range (δ_(ml)and δ_(mu)) for the least squares approximation and the color conversionmodel are optimized by letting δ_(ml)=0, δ_(mu)=0, and k=0, such thats_(k)={(L, a, b)_(jk)|7≦C*_(abjk)≦30}, and initiating the clusteringalgorithm to find the clusters. Then, R3→R is optimized for the CIELabto C, M, and Y conversions, respectively. The multidimensional Newton'smethod of R3→R3 is performed for the CIELab to CMY conversion and themean error and maximum error are recorded.

Afterwards, let δ_(ml)=2.5, δ_(mu)=2.5, and k=0, such that s_(k)=(L, a,b)_(jk)|7−2.5k≦C*_(abjk)≦30+2.5k,k=1,2, . . . , repeat the optimizationprocess stopping when all data points are included in s_(k). The leastsquare errors are obtained by calculating the differences between theoriginal values and the estimated values. The optimized models indicatethe Gaussian function for the cyan component, the inverse multiquadraticfor the magenta component, and the multiquadratic function for theyellow component, thereby providing inside gamut color conversions fromCIELab to CMY, g₀, and CMY to CIELab, f₀.

Overdetermined Linear Systems in the L_(∞) Sense

In saturated color areas, the perceptual color difference tolerance isgenerally higher than the magnitude of machine and instrument errors.Data outliers do not play a significant role in the color conversion.The maximum error is weighted in the color conversion qualityassessments. Thus, it is desirable to reduce the maximum approximationerror by solving the problem via the L_(∞) approximation, also referredas the Chebyshev norm. From norm inequalities if L_(∞) norm iscontrolled, then the L₁ and L₂ error is generally small.

Many algorithms have been developed for this type of applicationfunction approximation using the L_(∞) norm that are locally equivalentto the exchange algorithm. However, there is limited research onscattered data fitting techniques using the L_(∞) norm especially forthe large sparse problems. Large sparse problems are often converted tolinear programming problems using the linearized sub-problem as thebasis. Among these algorithms, the Barrodale and Phillips dual algorithm(the “BP” algorithm) and the Bartels, Conn, and Li primal algorithm (the“BC” algorithm) are suitable algorithms for solving overdeterminedlinear systems Rn→R in the Chebyshev norm.

The BP and BC algorithms may be suitable in view of the similarity ofthe formulations with the RBF approach. As the BP and BC algorithms aregenerally used in R^(n)→R, the extension of the BP dual algorithm may beextended to R3→R3 for color conversion from CIELab to CMY. In thisregard, the Chebyshev norm may be described by letting

Ax=b,   (eq. 41)

where A=[a₁, a₂, . . . , a_(n)] ε Rmxm, m>n≧2, and b^(T)=[β₁, β₂, . . ., β_(m)] ε Rm. Then, the vector x ε Rn is determined such that

$\begin{matrix}{{{{Ax} - b}}_{\infty} = {\max\limits_{1 \leq i \leq m}{{{a_{i}^{T}x} - \beta_{i}}}}} & \left( {{Eq}.\mspace{14mu} 42} \right)\end{matrix}$

is minimized. L is then defined as a linear space spanned by thefunctions {φ_(i)}

L=<φ ₁(z),φ₂(z), . . . , φ_(n)(z)>,   (Eq. 43)

where each φ_(i) is continuous on [a,b],j=1, . . . , n. A function φ ε Lmay then be expressed as

$\begin{matrix}{{{\varphi (z)} = {{\sum\limits_{j = 1}^{n}{\alpha_{j}{\varphi_{j}(z)}}} \in L}},} & \left( {{Eq}.\mspace{14mu} 44} \right)\end{matrix}$

where z ε [a, b]. The RBF is a relatively good example of such a space.Given m>n≧2, data points {(z_(i); y_(i)); i=1; 2; . . . ; m}, a=(a₁, . .. , a_(n)) is determined to minimize

$\begin{matrix}{\begin{matrix}{{e\left( {\alpha_{1},\ldots \mspace{11mu},\alpha_{n}} \right)} = {{y_{i} - {\sum\limits_{j = 1}^{n}{\alpha_{j}{\varphi_{j}\left( z_{i} \right)}}}}}_{\infty}} \\{= {\max\limits_{1 \leq i \leq m}{{{y_{i} - {\sum\limits_{j = 1}^{n}{\alpha_{j}{\varphi_{j}\left( z_{i} \right)}}}}}.}}}\end{matrix}{{Then},}} & \begin{matrix}\begin{matrix}\left( {{Eq}.\mspace{14mu} 45} \right) \\\;\end{matrix} \\\begin{matrix}\; \\\left( {{Eq}.\mspace{14mu} 46} \right)\end{matrix}\end{matrix} \\{x = {\alpha = \left\{ {a_{1},\alpha_{2},\ldots \mspace{11mu},\alpha_{n}} \right\}}} & \left( {{Eq}.\mspace{14mu} 47} \right) \\{{A = \begin{pmatrix}{\varphi_{1}\left( z_{1} \right)} & {\varphi_{2}\left( z_{1} \right)} & \ldots & {\varphi_{n}\left( z_{1} \right)} \\\ldots & \; & \; & \; \\\ldots & \; & \; & \; \\{\varphi_{1}\left( z_{m} \right)} & {\varphi_{2}\left( z_{m} \right)} & \ldots & {\varphi_{n}\left( z_{m} \right)}\end{pmatrix}}{and}} & \left( {{Eq}.\mspace{14mu} 48} \right) \\{b = {\left( {y_{1}\mspace{14mu} y_{2}\mspace{14mu} \ldots \mspace{14mu} y_{m}} \right)^{T}.}} & \left( {{Eq}.\mspace{14mu} 49} \right)\end{matrix}$

Next, the L_(∞) approximation may be converted to a standard form of thelinear programming. First, however, the primal form of the problem isestablished by letting ξ denote the maximum error. For example,

$\begin{matrix}{{\xi = {\max\limits_{1 \leq i \leq m}{e_{i}\left( {\alpha_{1},\alpha_{2},\ldots \mspace{11mu},\alpha_{n}} \right)}}},} & \left( {{Eq}.\mspace{14mu} 50} \right)\end{matrix}$

where e is defined in Eq. 41.

The primal problem may then be formulated by minimizing ξ subject to

$\begin{matrix}{= {{\xi + {\sum\limits_{j - 1}^{n}\; {\alpha_{j}{\varphi_{j}\left( z_{i} \right)}}}} \geq {y_{i}\mspace{14mu} {and}}}} & \left( {{Eq}.\mspace{14mu} 51} \right) \\{{= {{\xi - {\sum\limits_{j - 1}^{n}\; {\alpha_{j}{\varphi_{j}\left( z_{i} \right)}}}} \geq {- y_{i}}}},} & \left( {{Eq}.\mspace{14mu} 52} \right)\end{matrix}$

for every i where

ξ≧0   (Eq. 53)

and

a₁, a₂, . . . , a_(n) are unrestricted.   (Eq. 54)

The standard form of the corresponding dual linear program may betransformed from the above primal linear program via

$\begin{matrix}{\max {\sum\limits_{i = 1}^{m}\; {\left( {{\underset{\_}{\sigma}}_{i} - {\underset{\_}{\tau}}_{i}} \right)y_{i}}}} & \left( {{Eq}.\mspace{14mu} 55} \right)\end{matrix}$

subject to

$\begin{matrix}{{\sum\limits_{i = 1}^{m}\; {\left( {{\underset{\_}{\sigma}}_{i} - {\underset{\_}{\tau}}_{i}} \right)\varphi \; j_{i}}} = 0} & \left( {{Eq}.\mspace{14mu} 56} \right) \\{{\sum\limits_{i = 1}^{m}\; \left( {{\underset{\_}{\sigma}}_{i} + {\underset{\_}{\tau}}_{i}} \right)} \leq 1} & \left( {{Eq}.\mspace{14mu} 57} \right) \\{and} & \; \\{{\underset{\_}{\sigma}}_{i},{{\underset{\_}{\tau}}_{i} \geq 0},{i = 1},2,{\ldots \mspace{11mu} m},} & \left( {{Eq}.\mspace{14mu} 58} \right) \\{{{where}\mspace{14mu} \varphi_{ji}} = {{\varphi_{j}\left( z_{i} \right)}.}} & \;\end{matrix}$

To understand the BP dual algorithm, an understanding of the exchangealgorithm and related theorems which almost all Chebyshev approximationalgorithms are based is helpful. In this regard, A may be defined as alinear subspace of C [a,b], f ε C [a,b], and L is any closed subset of[a,b], p ε A, such that

e*(x)=f(x)−p*(x)   (Eq. 59)

and

L _(M) ={x ε L|e*(x) obtains extreme values}.   (eq. 60)

-   p* provides for a minimax solution to be approximation from A to f,    if ∀ p ε A, such that

e*(x)p(x)=[f(x)−p*(x)]p(x)>0, ∀x ε L _(M)   (Eq. 61)

(e.g., e*(x) and p(x) are of same signs at all extreme points).

In this regard, the “Haar” condition defines a class of functions uponwhich the Exchange Algorithm can be applied. For example, A may bedefined as an (n+1) dimensional linear subspace of C[a, b]. A satisfiesthe Haar condition if ∀ p ε A, the number of roots of the equation

{p(x)=0;a≦x≦b}  (Eq. 62)

is less than the dimension of A.

A derived condition states that if {ξ_(i):i=1; 2, . . . , k} ⊂ [a, b],k≦n then ∃ p ε A, such that p changes signs at ξ_(i):i=1; 2, . . . , kand p has no other roots. The RBFs in the form of Eq. 11 meet the Harrcondition because of the linear part in the equation. Then, let A be an(n+1)-dimensional linear subspace of C[a, b] that satisfies the Haarcondition and let f be any function in C[a, b]. Then p* is the minimaxapproximation from A to f, if there exists (n+2) points {ξ*_(i);i=0; 1;. . . ; n+1}, such that the conditions

a≦ξ* ₀ξ*₁< . . . <ξ*_(n+1) ≦b,   1.

|f(ξ*_(i))−p*(ξ*_(i))|=∥f−p*∥ _(∞) ,i=0,1, . . . , n+1 and   2.

f(ξ*_(i+1))−p*(ξ*_(i+1))=−[f(ξ*_(i))−p*(ξ*_(i))],i=0, . . . n,alternative signs]] are obtained.   3.

Then, let A be a finite-dimension subspace of C[a, b] that satisfies theHaar condition and let {x_(i):i=1; 2; . . . ;m} be a set of distinctpoints from [a, b], where m is not less than the dimension of A. For anyf in C[a, b], let the one-point exchange algorithm be applied tocalculate the element of A that minimizes expression

$\begin{matrix}{{\max\limits_{{i - 1},2,{\ldots \mspace{11mu} m}}{{{f\left( x_{i} \right)} - {p\left( x_{i} \right)}}}},{p \in A}} & \left( {{Eq}.\mspace{14mu} 63} \right)\end{matrix}$

Then, the required approximation to f is obtained in a finite number ofiterations.

In color conversion may compose a set of discrete data points. Thus, theExchange Algorithm may be based on the Minimax Characterization Theorem(MCT) as shown in FIG. 4. For each iteration k, the MCT adjusts areference {ξ_(i):i=0; 1; . . . ; n+1} and a trial approximation p_(k) isdetermined. On a discrete set of may be known on a set of points {x_(i);i=1; 2; : : : ; m}. To determine a trial approximation p_(k) in eachiteration, let

h=f(ξ₀)−p _(k)(ξ₀)   (Eq. 64)

and solve

f(ξ_(i))−p(ξ_(i))=(−1)^(i) h,i=0,1,2, . . . , n+1   (Eq. 65)

The iteration stops when p_(k) satisfies the conditions in the MCT.

For solving the dual problem, the BP dual algorithm is reviewed. The BPdual algorithm is based on the standard form of the simplex methodcomposed of three stages and introduces artificial variables a_(j) inEq. 55 and a slack variable ξ in Eq. 56. The above dual problem therebybecomes

$\begin{matrix}{\max {\sum\limits_{i = 1}^{m}\; {\left( {{\underset{\_}{\sigma}}_{i} - {\underset{\_}{\tau}}_{i}} \right)y_{i}}}} & \left( {{Eq}.\mspace{14mu} 66} \right)\end{matrix}$

subject to

$\begin{matrix}{{{{\sum\limits_{i = 1}^{m}\; {\left( {{\underset{\_}{\sigma}}_{i} - {\underset{\_}{\tau}}_{i}} \right)\varphi_{j_{i}}}} + {\underset{\_}{\alpha}}_{i}} = 0},{j = 1},2,\ldots,n} & \left( {{Eq}.\mspace{14mu} 67} \right) \\{{{\sum\limits_{i = 1}^{m}\; \left( {{\underset{\_}{\sigma}}_{i} + {\underset{\_}{\tau}}_{i}} \right)} + \underset{\_}{\xi}} = 1} & \left( {{Eq}.\mspace{14mu} 68} \right) \\{\xi,{\underset{\_}{\alpha}}_{i},{{- {\underset{\_}{\tau}}_{i}} \geq 0},{j = 1},2,\ldots \mspace{11mu},{m.}} & \left( {{Eq}.\mspace{14mu} 69} \right)\end{matrix}$

Denoting the column vectors formed by the coefficients of σ _(i), α_(i), τ _(i), and ξ by s _(i), t _(i), α _(j), and w, then

si+t _(i)=2w   (Eq. 70)

Then, let b denote the right-hand side of Eq. 52 such that

b=w=(0 0 . . . 0 1)^(T) =e ^(T) _(n+1).   (Eq. 71)

The first two stages of the BP dual algorithm include determining anefficient starting point for the third stage, which is equivalent to theexchange algorithm for the linear minimax approximation. In doing so,the BP dual algorithm uses the table below:

Then, Φ is denoted as an m×n matrix:

$\begin{matrix}{\Phi = {\begin{pmatrix}{\varphi_{1,1}\varphi_{2,1}\ldots \mspace{11mu} \varphi_{n,1}} \\{\varphi_{1,2}\varphi_{2,2}\ldots \mspace{11mu} \varphi_{n,2}} \\\ldots \\\ldots \\{\varphi_{1,m}\varphi_{2,m}\ldots \mspace{11mu} \varphi_{n,m}}\end{pmatrix}\mspace{14mu} {and}}} & \left( {{Eq}.\mspace{14mu} 72} \right)\end{matrix}$

the rank of Φ is denoted as k.

In the first stage, the first k simplex iterations are performed. Thealgorithm finds the maximum residual at each iteration and reduce it tozero. The nonbasic variables that can enter the basis to become basicvariables are σ _(i)'s, and basic variables that can leave the basis tobecome nonbasic variables are restricted to α _(i)s. Basic variable ξremains zero in the basis. At the end of stage 1, there are k residuesat zero. σ _(t) that enters the basis corresponds to that with thelargest absolute reduced cost (e.g., choosing the largest residual inthe primal problem). The pivot may be selected from this column as thelargest absolute value corresponding to an α _(j) in the basis. At thet^(th) iteration, with the corresponding σ _(t) entering the basis,there are t zero residuals. For example,

$\begin{matrix}{{y_{t} - {\sum\limits_{j = 1}^{n}\; {\alpha_{j}{\varphi_{j}\left( z_{t} \right)}}}} = {\xi = 0.}} & \left( {{Eq}.\mspace{14mu} 73} \right)\end{matrix}$

In the second stage, the (k+1)^(th) simplex iteration is performed andthe variable is forced to leave the basis. This means the error ξ in theprimal problem is increased above zero. At the end of the second stage,the k+1 residuals of magnitude equal to resulting value of ξ. The(k+1)^(th) pivotal column may be chosen also corresponding to thelargest absolute marginal cost (e.g., equivalent to choosing the largestresidual in the primal problem). The σ_(j)'s and τ_(j)'s are exchangedif the marginal cost is negative. If there are positive values in thefirst n elements in the pivotal column, these values are changed tonegative by adding twice of the corresponding rows to the (n+1)^(th)pivotal row and changing the sign of those original rows.

The third stage is equivalent to the exchange algorithm for a linearminimax approximation. The pivotal column may be chosen corresponding tothe most negative marginal cost or a ratio selection rule. The iterationcontinues until all marginal cost is nonnegative. Every iteration mayincrease the value of ξ until ξ=e*, where e* is the minimum value in Eq.45.

Alternatively, the BC algorithm may be used to determine the startingpoint for color conversion. In doing so, the primal formulation in Eq.45 is rewritten by letting

$\begin{matrix}{c_{0} = \begin{pmatrix}1 & 0 & 0 & 0 & \ldots & 0\end{pmatrix}_{n + 1}^{T}} & \left( {{Eq}.\mspace{14mu} 74} \right) \\{c_{i} = \begin{pmatrix}1 \\{- a_{i}}\end{pmatrix}_{n + 1}} & \left( {{Eq}.\mspace{14mu} 75} \right) \\{c_{m - i} = \begin{pmatrix}1 \\{+ a_{i}}\end{pmatrix}_{n + 1}} & \left( {{Eq}.\mspace{14mu} 76} \right) \\{\delta_{i} = {- y_{i}}} & \left( {{Eq}.\mspace{14mu} 77} \right) \\{{\delta_{m + i} = y_{i}},{i = 1},2,\ldots,m,{and}} & \left( {{Eq}.\mspace{14mu} 78} \right) \\{v = \begin{pmatrix}\xi \\\alpha\end{pmatrix}_{n + 1}} & \left( {{Eq}.\mspace{14mu} 79} \right)\end{matrix}$

The above primal problem is then converted to an unconstrainedoptimization problem via a piecewise linear penalty function. Forexample, let μ>0 be a fixed parameter and define

$\begin{matrix}{{p\left( {{vi},\mu} \right)} = {{\mu \; c_{0}^{T}v} - {\sum\limits_{j = 1}^{2m}\; {\min \left( {0,{{c_{j}^{T}v} - \delta_{j}}} \right)}}}} & \left( {{Eq}.\mspace{14mu} 80} \right)\end{matrix}$

For any arbitrary v ε Rn+1, Eq. 80 can be expanded into

$\begin{matrix}{{p\left( {v,\mu} \right)} = {{\mu \; c_{0}^{T}v} - {\sum\limits_{j \in I^{0}}\; {\min \left( {0,{{c_{j}^{T}v} - \delta_{j}}} \right)}}}} & \left( {{Eq}.\mspace{14mu} 81} \right) \\{- {\sum\limits_{j \in I^{+}}\; {\min \left( {0,{{c_{j}^{T}v} - \delta_{j}}} \right)}}} & \left( {{Eq}.\mspace{14mu} 82} \right) \\{{- {\sum\limits_{j \in I^{-}}\; {\min \left( {0,{{c_{j}^{T}v} - \delta_{j}}} \right)}}},} & \left( {{Eq}.\mspace{14mu} 83} \right) \\{where} & \; \\\begin{matrix}{{I^{0} = {{I^{0}(v)} = {\left\{ {{j{c_{j}^{T}v}} = \delta_{j}} \right\} = \left\{ {j_{i},j_{2},\ldots,j_{k}} \right\}}}},} \\{{I^{+} = {{I^{+}(v)} = \left\{ {j{{c_{j}^{T}v} > \delta_{i}}} \right\}}},{and}} \\{I^{-} = {{I^{-}(v)} = {\left\{ {j{{c_{j}^{T}v} < \delta_{j}}} \right\}.}}}\end{matrix} & \;\end{matrix}$

For any d ε

and any λ≦0 sufficiently small,

$\begin{matrix}{{p\left( {{v + {\lambda \; d}},\mu} \right)} = {{p\left( {v + \mu} \right)} + {\lambda\left\lbrack {{\mu \; c_{0}^{T\;}d} - {\sum\limits_{j \in I^{-}}^{\;}\; {c_{j}^{T}d}} - {\sum\limits_{j \in I^{0}}^{\;}\; {\min \left( {0,{c_{j}^{T}d}} \right)}}} \right\rbrack}}} & \left( {{Eq}.\mspace{14mu} 84} \right) \\{\mspace{130mu} {{= {{p\left( {v + \mu} \right)} + {\lambda \; h^{T\;}d} + {\lambda {\sum\limits_{j \in I^{0}}^{\;}\; {\sigma_{j}^{-}c_{j}^{T}d}}}}},}} & \left( {{Eq}.\mspace{14mu} 85} \right) \\{where} & \; \\{{h = {{\mu \; c_{0}} - {\sum\limits_{j \in I^{-}}^{\;}\; c_{j}}}},} & \; \\{and} & \; \\{\sigma_{j}^{-} = {{\sigma_{j}^{-}(d)} = \left\{ {{\begin{matrix}0 & {{{if}\mspace{14mu} c_{j}^{T}d} \geq 0} \\1 & {{{{if}\mspace{14mu} c_{j}^{T}d} < 0},}\end{matrix}j} \in {I^{0}.}} \right.}} & \;\end{matrix}$

Then, define the matrix

N=[c _(j1) , c _(j2) , . . . c _(jk)]  (Eq. 86)

and let

P=I−N(N ^(T) N)⁻¹ N ^(T).   (Eq. 87)

The matrix P is the orthogonal projector onto the null space of N^(T).To define search direction d such that p(v+λd, μ)<p(v, μ),

d=−Ph   (Eq. 88)

is established. When Ph≠0, apply d defined in Eq. 88 to Eq. 85 such that

$\begin{matrix}\begin{matrix}{{c_{j}^{T}d} = {{- c_{j}^{T}}{Ph}}} \\{= {{- c_{j}^{T}}{PPh}}} \\{= {{- \left( {Pc}_{j} \right)^{T}}{Ph}}} \\{= 0}\end{matrix} & \left( {{Eq}.\mspace{14mu} 89} \right)\end{matrix}$

Since P is the orthogonal projector onto the null space of N^(T). TheEq. 84 becomes:

p(v+λd, μ)=p(v, μ)+λh ^(T) d   (Eq. 90)

h^(T)d<0, and   (Eq. 91)

p(v+λd, μ)<p(v, μ).   (Eq. 92)

Thus, d serves as a descent direction for the function p.

When Ph=0, assuming columns of N are linearly independent, then

$\begin{matrix}{h = {\sum\limits_{i = 1}^{k}\; {\eta_{i}{c_{ji}.}}}} & \left( {{Eq}.\mspace{14mu} 93} \right)\end{matrix}$

and Eq. 84 becomes

$\begin{matrix}{{p\left( {{v + {\lambda \; d}},\mu} \right)} = {{p\left( {v,\mu} \right)} + {\lambda {\sum\limits_{i = 1}^{k}\; {\left\lbrack {\eta_{i} + \sigma_{ji}^{-}} \right\rbrack c_{j}^{T}{d.}}}}}} & \left( {{Eq}.\mspace{14mu} 94} \right)\end{matrix}$

In the sub case of a: η_(i*)<0, for some i* ε I⁰, d is chosen such that

c^(T) _(j)d=0, j ε I⁰, j≠i,   (Eq. 95)

c_(i) ^(T)*d=1   (Eq. 96)

This implies that σ_(j) ⁻=0<0, j ε I⁰,. Because

h ^(T) d=η _(i*) c _(i*) ^(T) d=η<0,   (Eq. 97)

then

p(v+λd, μ)<p(v, μ).   (Eq. 98)

In the sub case of b: η_(i*)<0 for all i* εI ⁰, if no constraints areviolated, it is the optimal. Otherwise, reduce μ and start anotheriteration. For example, {choose step size λ}

$\begin{matrix}{{\Lambda^{{(1)}^{\underset{\_}{\underset{\_}{def}}}}\left\{ {{{\lambda \lambda} = \frac{\delta_{j} - {c_{j}^{T}v}}{c_{j}^{T}d}},{\lambda > 0},{j \in {I^{+}\bigcup I^{-}}}} \right\}},} & \left( {{Eq}.\mspace{14mu} 99} \right) \\{\lambda = {\min \; \Lambda^{(1)}}} & \left( {{Eq}.\mspace{14mu} 100} \right)\end{matrix}$

and update v=v+λd.

As an example, the L_(∞) norm approximation may be applied to colorconversions from CMYK to CIELab and CIELab to CMY. The formulation ofthis color conversion in the primal form is described in the followingusing the CIELab to CMY color conversion by minimizing

ξ=ξ_(c)+ξ_(m)ξ_(y)   (Eq. 101)

subject to

$\begin{matrix}{{\xi_{c} + {\sum\limits_{j = 1}^{n + 4}\; {\alpha_{cj}{F_{cj}\left( z_{i} \right)}}}} \geq c_{i}} & \left( {{Eq}.\mspace{14mu} 102} \right) \\{{\xi_{c} - {\sum\limits_{j = 1}^{n + 4}\; {\alpha_{cj}{F_{cj}\left( z_{i} \right)}}}} \geq {- c_{i}}} & \left( {{Eq}.\mspace{14mu} 103} \right) \\{{\xi_{m} + {\sum\limits_{j = 1}^{n + 4}\; {\alpha_{mj}{F_{mj}\left( z_{i} \right)}}}} \geq m_{i}} & \left( {{Eq}.\mspace{14mu} 104} \right) \\{{\xi_{m} - {\sum\limits_{j = 1}^{n + 4}\; {\alpha_{mj}{F_{mj}\left( z_{i} \right)}}}} \geq {- m_{i}}} & \left( {{Eq}.\mspace{14mu} 105} \right) \\{{\xi_{y} + {\sum\limits_{j = 1}^{n + 4}\; {\alpha_{yj}{F_{yj}\left( z_{i} \right)}}}} \geq y_{i}} & \left( {{Eq}.\mspace{14mu} 106} \right) \\{{{\xi_{y} - {\sum\limits_{j = 1}^{n + 4}\; {\alpha_{yj}{F_{yj}\left( z_{i} \right)}}}} \geq {- y_{i}}},} & \left( {{Eq}.\mspace{14mu} 107} \right) \\{where} & \; \\{{\xi \geq 0},} & \;\end{matrix}$

-   a₁, a₂, . . . a_(n+4) are unrestricted,-   z_(i)={L, a, b}i, i=1, . . . , m,-   F_(ji)=[1, z_(i), φ₁(z_(i)), φ₂(z_(i)), . . . , φ_(n)(z_(i))]_(4+n),-   φ_(ji)=φ_(j)(z_(i)),-   φ_(j) is one of the RBFs, n is the number of clusters, and m is the    number of color data points.

The problem may then be rewritten by minimizing c₀ ^(T)v subject to

$\begin{matrix}{{{c_{j}^{T}v} \geq \delta_{j}},{j = 1},2,{\ldots \mspace{11mu} 6m}} & \left( {{Eq}.\mspace{14mu} 108} \right) \\{{where}{c_{0} = \begin{matrix}\left( 1 \right. & 0 & \ldots & 1 & 0 & \ldots & 1 & 0 & \ldots & \left. 0 \right)_{{3n} + 15}^{T}\end{matrix}}{v = \begin{pmatrix}\xi_{c} & \alpha_{c\; 1} & \ldots & \alpha_{cn} & \xi_{c} & \alpha_{m\; 1} & \ldots & \alpha_{mn} & \xi_{y} & \alpha_{y\; 1} \\\ldots & \alpha_{yn} & \; & \; & \; & \; & \; & \; & \; & \;\end{pmatrix}_{{3n} + 15}^{T}}{c_{i} = \begin{matrix}\left( 1 \right. & {- \alpha_{ci}} & \ldots & 0 & \ldots & \left. 0_{{2n} + 8} \right)_{{3n} + 15}^{T}\end{matrix}}{c_{m + i} = \begin{matrix}\left( 1 \right. & \alpha_{ci} & \ldots & 0 & \ldots & \left. 0_{{2n} + 8} \right)_{{3n} + 15}^{T}\end{matrix}}{c_{{2m} + i} = \begin{matrix}\left( 0 \right. & \ldots & 0_{n + 4} & \ldots & 1 & \ldots & {- \alpha_{mi}} & 0 & \ldots & \left. 0_{n + 4} \right)_{{3n} + 15}^{T}\end{matrix}}{c_{{3m} + i} = \begin{matrix}\left( 0 \right. & \ldots & 0_{n + 4} & \ldots & 1 & \alpha_{mi} & 0 & \ldots & \left. 0_{n + 4} \right)_{{3n} + 15}^{T}\end{matrix}}{c_{{4m} + i} = \begin{matrix}\left( 0 \right. & \ldots & 0_{{2n} + 8} & 1 & \left. {- \alpha_{yi}^{\;}} \right)_{{3n} + 15}^{T}\end{matrix}}{c_{{5m} + i} = \begin{matrix}\left( 0 \right. & \ldots & 0_{{2n} + 8} & 1 & \left. \alpha_{yi}^{\;} \right)_{{3n} + 15}^{T}\end{matrix}}{\delta_{i} = {- C_{i}}}{\delta_{m + i} = C_{i}}{\delta_{{2m} + i} = {- M_{i}}}{\delta_{{3m} + i} = M_{i}}{\delta_{{4m} + i} = {- Y_{i}}}{{\delta_{{5m} + i} = Y_{i}},{and}}{{i = 1},\ldots,{m.}}} & \;\end{matrix}$

The primal problem may then be divided into steps using the modified BCalgorithm. For example, the primal problem with constraints may beconverted to an unconstrained problem via a piecewise linear penaltyfunction by letting μ>0 be a fixed parameter. Next,

$\begin{matrix}{{p\left( {v,\mu} \right)} = {{\mu \; c_{0}^{T}v} - {\sum\limits_{j = 1}^{6m}\; {\min \left( {0,{{c_{j}^{T}v} - \delta_{j}}} \right)}}}} & \left( {{Eq}.\mspace{14mu} 109} \right)\end{matrix}$

is defined.For any arbitrary v ε

$\begin{matrix}{{{{p\left( {v,\mu} \right)} = {{\mu \; c_{0}^{T}v} - {\sum\limits_{j \in I^{0}}^{\;}\; {\min \left( {0,{{c_{j}^{T}v} - \delta_{j}}} \right)}} - {\sum\limits_{j \in I^{+}}\; {\min \left( {0,{{c_{j}^{T}v} - \delta_{j}}} \right)}} - {\sum\limits_{j \in I^{-}}\; {\min \left( {0,{{c_{j}^{T}v} - \delta_{j}}} \right)}}}},{where}}{{I^{0} = {{I^{0}(v)} = {\left\{ {{j{c_{j}^{T}v}} = \delta_{j}} \right\} = \left\{ {j_{1},j_{2},\ldots \;,\; j_{k}} \right\}}}},{I^{+} = {{I^{+}(v)} = \left\{ {j{{c_{j}^{T}v} > \delta_{j}}} \right\}}},{and}}{I^{-} = {{I^{-}(v)} = {\left\{ {j{{c_{j}^{T}v} < \delta_{j}}} \right\}.}}}} & \left( {{Eq}.\mspace{14mu} 110} \right)\end{matrix}$

Then, the sets I⁰, I⁺, and I⁻ are determined and the search directionand step size are updated.

As an example of color conversion in dual for CIELab to CMY,

$\begin{matrix}{{{Maximize}{\sum\limits_{i = 1}^{m}\; {\left( {{\underset{\_}{\sigma}}_{ci} - {\underset{\_}{\tau}}_{ci}} \right)c_{i}}}} + {\left( {{\underset{\_}{\sigma}}_{mi} - {\underset{\_}{\tau}}_{mi}} \right)m_{i}} + {\left( {{\underset{\_}{\sigma}}_{yi} - {\underset{\_}{\tau}}_{yi}} \right)y_{i}}} & \left( {{Eq}.\mspace{14mu} 111} \right)\end{matrix}$

subject to

$\begin{matrix}{{\sum\limits_{i = 1}^{m}\; {\left( {{\underset{\_}{\sigma}}_{ci} + {\underset{\_}{\tau}}_{ci}} \right)F_{C_{ji}}^{T}}} = 0} & \left( {{Eq}.\mspace{14mu} 112} \right) \\{{\sum\limits_{i = 1}^{m}\; \left( {{\underset{\_}{\sigma}}_{ci} - {\underset{\_}{\tau}}_{ci}} \right)} \leq 1} & \left( {{Eq}.\mspace{14mu} 113} \right) \\{{\sum\limits_{i = 1}^{m}\; {\left( {{\underset{\_}{\sigma}}_{Mi} + {\underset{\_}{\tau}}_{Mi}} \right)F_{M_{ji}}^{T}}} = 0} & \left( {{Eq}.\mspace{14mu} 114} \right) \\{{\sum\limits_{i = 1}^{m}\; \left( {{\underset{\_}{\sigma}}_{Mi} - {\underset{\_}{\tau}}_{Mi}} \right)} \leq 1} & \left( {{Eq}.\mspace{14mu} 115} \right) \\{{\sum\limits_{i = 1}^{m}\; {\left( {{\underset{\_}{\sigma}}_{Yi} - {\underset{\_}{\tau}}_{Yi}} \right)F_{Y_{ji}}^{T}}} = 0} & \left( {{Eq}.\mspace{14mu} 116} \right) \\{{\sum\limits_{i = 1}^{m}\; \left( {{\underset{\_}{\sigma}}_{Yi} + {\underset{\_}{\tau}}_{Yi}} \right)} \leq 1} & \left( {{Eq}.\mspace{14mu} 117} \right) \\{and} & \; \\{{\underset{\_}{\sigma}}_{i},{{\underset{\_}{\tau}}_{i} \geq 0},{i = 1},2,\ldots \;,m,} & \left( {{Eq}.\mspace{14mu} 118} \right) \\{{{where}{z_{i} = \left\{ {L,a,b} \right\}_{i}},{i = 1},2,\ldots \;,m}{F_{ji} = \left\lbrack {1,z_{i},{\phi_{1}\left( z_{i} \right)},{\phi_{2}\left( z_{i} \right)},\ldots \;,{\phi_{n}\left( z_{i} \right)}} \right\rbrack_{n + 4}}{{\phi_{ji} = {\phi_{j}\left( z_{j} \right)}},}} & \;\end{matrix}$

φ_(j) is one of the RBFs, n is the number of clusters, and m is thenumber of color data points.

The CIELab to CMY conversion may then be represented by the initialsimplex table:

and the rank of the matrix F_(n×m), k is thereby determined.

Thereafter, the first step BP simplex algorithm may be performed via ksimplex iterations for each of F_(Cn×m) for the cyan component, F_(Mn×m)for the magenta component, and F_(Yn×m) for the yellow component. Thevectors s_(C) are allowed to enter the basis and the vectors a_(C) canleave the basis for F_(Cn×m). The vectors s_(M) are allowed to enter thebasis and the vectors a_(M) can leave the basis for F_(Mn×m). Thevectors s_(Y) are allowed to enter the basis and the vectors a_(Y) canleave the basis for Y_(n×m).

The second step of the BP simplex algorithm may be performed by forcingξ_(C), ξ_(M), and ξ_(Y) to leave the basis. The corresponding vectorsentering the basis are s_(C), s_(M), and s_(Y) respectively. If thecorresponding marginal cost of s_(C), s_(M), or s_(Y) is positive, thensi is changed to t_(i) using

t _(i)=2w3_(n+4) −s _(i),

where

w _(C)=(0 0 . . . 0 1 0 0)^(T) =e _(3n+1) ^(T)

w _(M)=(0 0 . . . 0 0 1 0 0)^(T) =e _(3n|2) ^(T), and

w _(Y)=(0 0 . . . 0 0 0 1 0)^(T) =e _(3n+3) ^(T).

If s_(i) is the pivotal column, r is the row vector of pivotal row. Forany j^(th) row vector r_(j) containing s_(ij)>0, where j≠3n+1 forF_(Cn×m), j≠3n+2 for F_(Mn×m), and j≠3n+3 for F_(Yn×m), change r tor+2r_(j) and change r_(j) to −r_(j). The third step of the BP simplexiterations may then be performed until all the marginal costs arenon-negative.

Color Conversion using the L₁ Norm Approximation

Color conversion in neutral color areas plays an important role in colorreproduction quality. The human visual system (HVS) is sensitive tocolor changes in neutral color areas. Additionally, the color conversionerrors due to machine stability, paper uniformity, and instrumentaccuracy are more profound relative to the human color differencetolerances. Thus, there are more chances for the existences of outliersusing color conversion. A color conversion should minimize the negativeimpact of outliers. The L₁ norm approximation may assist in the regiondue to the robustness to outliers.

Among all the L₁ algorithms, the Barrodale and Roberts' primal algorithm(BR algorithm) is relatively robust. Similar to the problem setup forthe L_(∞), the overdetermined linear system of the general L₁approximation problem may be stated as follows:

Ax=b,   (Eq. 119)

where

A=[a ₁ , a ₂ , . . . , a _(n)] ε

(m>n≧2)

and

b^(T)=[β₁, β₂, . . . , β_(m)] ε

The L₁ norm approximation is used to determine vector x ε

such that

$\begin{matrix}{{{{Ax} - b}}_{1} = {\sum\limits_{i = 1}^{m}\; {\sum\limits_{j = 1}^{n}\; {{{a_{j}^{T}x_{i}} - \beta_{i}}}}}} & \left( {{Eq}.\mspace{14mu} 120} \right)\end{matrix}$

is minimized. In this regard, L is a linear space spanned by functionφ_(i)

L=<φ ₁(z), φ₂(z), . . . , φ_(n)(z)>,   (Eq. 121)

where each φ_(i) is continuous on [a, b], i=1, . . . , n. The function(P may be expressed as

$\begin{matrix}{{{\varphi (z)} = {{\sum\limits_{i = 1}^{n}\; {\alpha_{i}{\varphi_{i}(z)}}} \in L}},} & \left( {{Eq}.\mspace{14mu} 122} \right)\end{matrix}$

where z ε [a,b] ⊂

Given m>n≧2, data points (z_(i), y_(i)), a=(a₁, . . . a_(n)) isdetermined to minimize

$\begin{matrix}{{\left( {\alpha_{1}\ldots \mspace{11mu} \alpha_{n}} \right)} = {{y_{i} - {\sum\limits_{j = 1}^{n}\; {\alpha_{j}{\varphi_{j}\left( z_{i} \right)}}}}}_{1}} & \left( {{Eq}.\mspace{14mu} 123} \right) \\{\mspace{121mu} {= {\sum\limits_{j = 1}^{m}\; {{y_{j} - {\sum\limits_{j = 1}^{n}\; {\alpha_{j}{\varphi_{j}\left( z_{i} \right)}}}}}}}} & \left( {{Eq}.\mspace{14mu} 124} \right)\end{matrix}$

The equation then becomes:

$\begin{matrix}{x = a} & \; \\{A = \left( \begin{matrix}{\varphi_{1}\left( z_{1} \right)} & {\varphi_{2}\left( z_{1} \right)} & \ldots & {\varphi_{n}\left( z_{1} \right)} \\\ldots & \; & \; & \; \\\ldots & \; & \; & \; \\{\varphi_{1}\left( z_{m} \right)} & {\varphi_{2}\left( z_{m} \right)} & \ldots & {\varphi_{n}\left( z_{m} \right)}\end{matrix} \right.} & \left( {{Eq}.\mspace{14mu} 125} \right) \\{b = \left( {y_{1}y_{2}\ldots \mspace{11mu} y_{m}} \right)^{T}} & \;\end{matrix}$

In the RBF form, the function becomes

$\begin{matrix}{{{\varphi (z)} = {\alpha_{0} + {\sum\limits_{p = 1}^{k}\; {\alpha_{p}{z(p)}}} - {\sum\limits_{j = {n + k}}^{n + k}\; {\alpha_{j}{\varphi \left( {{z - c_{j}}} \right)}}}}},} & \left( {{Eq}.\mspace{14mu} 126} \right) \\{{{where}\mspace{14mu} x}\; = {\left( {\alpha_{0},\alpha_{1},\ldots \;,\alpha_{k},\alpha_{k + 1},\ldots \;,\alpha_{n + k}} \right)\mspace{14mu} {and}}} & \; \\{A = {\begin{pmatrix}1 & {z(1)} & {z_{1}(2)} & \ldots & {z_{1}(k)} & {\varphi_{1}\left( z_{1} \right)} & {\varphi_{2}\left( z_{1} \right)} & \ldots & {\varphi_{n}\left( z_{1} \right)} \\\ldots & \; & \; & \; & \; & \; & \; & \; & \; \\\ldots & \; & \; & \; & \; & \; & \; & \; & \; \\\ldots & \; & \; & \; & \; & \; & \; & \; & \; \\1 & {z_{m}(1)} & {z_{m}(2)} & \ldots & {z_{m}(k)} & {\varphi_{1}\left( z_{m} \right)} & {\varphi_{2}\left( z_{m} \right)} & \ldots & {\varphi_{n}\left( z_{m} \right)}\end{pmatrix}.}} & \;\end{matrix}$

The L₁ and L_(∞) norms in C[a, b] and in R_(n) are not strictly convex.Accordingly, A is a linear subspace of either the L₁ or the L_(∞) normedlinear space, as the uniqueness of best approximations from A to fdepends on properties of A and f. p* in A may be the better L₁approximation function to f A sign function is thus defined as:

$\begin{matrix}{{s*(x)} = \left\{ \begin{matrix}{{- 1},} & {{f(x)} < {p^{*}(x)}} \\{0,} & {{f(x)} = {p^{*}(x)}} \\{1,} & {{f(x)} > {p^{*}(x)}}\end{matrix} \right.} & \left( {{Eq}.\mspace{14mu} 127} \right)\end{matrix}$

A sufficient condition for the function p* is the L₁ approximation fromA to f which may be defined by letting A be a linear subspace of C[a,b], f be any function in C[a, b], and p* be any element of A, such thatthe set

L={x:f(x)=p*(x), a≦x≦b}  (Eq. 128)

is either empty or is composed of a finite number of intervals anddiscrete points. Then p* provides an L₁ norm approximation form A to f,if the inequality

$\begin{matrix}{{\left. {\int_{a}^{b}{{s^{*}(x)}{p(x)}\ {x}}} \right\} } \leq {\int_{L}^{\;}{{{p(x)}}\ {x}\mspace{20mu} (4.7)}}} & \left( {{Eq}.\mspace{14mu} 129} \right)\end{matrix}$

is satisfied for all p in A, where s* is the function of Eq. 127.

For discrete L₁ approximation, let the function values {f(x); t=1, 2, .. . , m} and the fixed positive weights be {w_(t); t=1, 2, . . . , m}. Ais a linear space of functions that may be defined on the point set{x_(t); t=1, 2, . . . , m}, and p* is any element of A. Then, ε containsthe points of {x_(t); t=1, 2, . . . , m} that satisfy the condition

p*(x _(t))=f(x _(t))   (Eq. 130)

and s* is the sign function

$\begin{matrix}{{s^{*}(x)} = \left\{ \begin{matrix}{{- 1},} & {{f\left( x_{t} \right)} < {p^{*}\left( x_{t} \right)}} \\{0,} & {{f\left( x_{t} \right)} = {p^{*}\left( x_{t} \right)}} \\{1,} & {{f\left( x_{t} \right)} > {p^{*}\left( x_{t} \right)}}\end{matrix} \right.} & \left( {{Eq}.\mspace{14mu} 131} \right)\end{matrix}$

t=1, 2, . . . , m. p* becomes the function in A that minimizes theexpression

$\begin{matrix}{{\sum\limits_{t = 1}^{m}\; {w_{t}{{{f\left( x_{t} \right)} - {p\left( x_{t} \right)}}}}},{p \in A}} & \left( {{Eq}.\mspace{14mu} 132} \right)\end{matrix}$

if the inequality

$\begin{matrix}{{{\sum\limits_{t = 1}^{m}\; {w_{t}{s^{*}\left( x_{t} \right)}{p\left( x_{t} \right)}}}} \leq {\sum\limits_{x_{t \in L}}^{m}\; {w_{t}{{p\left( x_{t} \right)}}}}} & \left( {{Eq}.\mspace{14mu} 133} \right)\end{matrix}$

holds for all p in A.

In order to improve the discrete L₁ approximation p*, the inequality inEq. 133 is tested to be satisfied for every element p in A. An optimalfunction p* can be obtained by searching for suitable interpolationpoints in L. For example, let the function values {f(x_(t)); t=1, 2, . .. , m} and the fixed positive weights be {w_(t); t=1, 2, . . . , m}.Then, A is a linear subspace of R^(m), where the component of eachvector p in A have the values {p(x_(t)); t=1, 2, . . . , m}. An elementp* in A minimizes Eq. 132, and the property that the zero vector is theelement p in A that satisfies the conditions {p(x_(t))=0; x_(t)*L},where the set L is defined in Eq. 128

With this is mind, the discrete L₁ approximation calculation becomes alinear programming problem. For example with φ={φ_(i); i=0, 1, . . . ,n} as a basis of the space A of approximation, the Eq. 133 becomes

$\begin{matrix}{\sum\limits_{t = 1}^{m}\; {w_{t}{{y_{t} - {\sum\limits_{i = 0}^{n}\; {\lambda_{i}{\varphi_{i}\left( x_{t} \right)}}}}}}} & \left( {{Eq}.\mspace{14mu} 134} \right)\end{matrix}$

where y_(t) is the output value corresponding to each x_(t). Define thebounds {μ_(t)≧0} and {u_(t)≧0} may be defined such that

$\begin{matrix}{{- v_{t}} \leq {y_{t} - {\sum\limits_{i = 0}^{n}\; {\lambda_{i}{\varphi_{i}\left( x_{t} \right)}}}} \leq u_{t}} & \left( {{Eq}.\mspace{14mu} 135} \right)\end{matrix}$

for {t=1, 2, . . . , m}. Then the linear programming problem is

$\begin{matrix}{{Minimize}{\sum\limits_{t = 1}^{m}\; {w_{t}\left( {u_{t} + v_{t}} \right)}}} & \left( {{Eq}.\mspace{14mu} 136} \right)\end{matrix}$

subject to

$\begin{matrix}{{{- v_{t}} \leq {y_{t} - {\sum\limits_{i = 0}^{n}\; {\lambda_{i}{\varphi_{i}\left( x_{t} \right)}}}} \leq u_{t}},{t = 1},2,{\ldots \mspace{11mu} m}} & \left( {{Eq}.\mspace{14mu} 137} \right) \\{{{ut} \geq 0}{{vt} \geq 0}{{{{for}\mspace{14mu} t} = 1},2,\ldots \;,{m.}}} & \;\end{matrix}$

The algorithm may be described as:

$\begin{matrix}{{{f_{i} - {\sum\limits_{j = 1}^{n}\; {\alpha_{j}\varphi_{j,i}}}} = {u_{i} - v_{i}}},} & \left( {{Eq}.\mspace{14mu} 138} \right)\end{matrix}$

where {u_(i), v_(i): i=1, 2, . . . , m} are nonnegative variables. Thenonnegative variables {b_(j), c_(j): j=1, 2, . . . , n} and the weightsa_(j)=b_(j)−c_(j) may then be defined. The primal problem becomes

$\begin{matrix}{{Minimizing}{\sum\limits_{i = 1}^{m}\; \left( {u_{i} + v_{i}} \right)}} & \left( {{Eq}.\mspace{14mu} 139} \right)\end{matrix}$

subject to

$\begin{matrix}{{f_{i} = {{\sum\limits_{j = 1}^{n}\; {\left( {b_{j} - c_{j}} \right)\varphi_{j,i}}} + u_{i} - v_{i}}},{i = 1},2,\ldots \;,m} & \left( {{Eq}.\mspace{14mu} 140} \right) \\{and} & \; \\{b_{j},c_{j},u_{i},{v_{i} \geq 0.}} & \left( {{Eq}.\mspace{14mu} 141} \right)\end{matrix}$

The algorithm may be implemented in two stages: For example, for thefirst n iterations, the pivotal columns are restricted to b_(j) andc_(j). The vector entering the basis is the one with the largestnonnegative marginal cost,

$\left( {{e.g.},{\sum\limits_{i = 1}^{m}\; \varphi_{j,i}}} \right).$

the vector leaving the basis is chosen among v_(i) and u_(i) for the onecausing the maximum reduction in the objective function. Next, thenonbasic u_(i) or v_(i) are interchanged with the basic u_(i) or v_(i).Neither b_(j) or c_(j) is allowed to leave the basis. The criteria forselecting vectors leaving and entering the basis are the same as that inthe first stage. Otherwise, the basic vectors b_(j) or c_(j) areinterchanged with the corresponding nonbasic vectors c_(j) or b_(j). Onmaking this algorithm efficient, the vector chosen to leave the basis isthe one that causes the maximum reduction in the objective function.

On the Multidimensional Neutral Color Conversion Problem, the BRalgorithm may be applied to the neutral color conversion by lettingz_(i)={(L, a, b)_(i)i=1, 2, . . . , m}. As an example, the colorconversion from CIELab to CMY is formulated:

$\begin{matrix}{{{Minimizing}\mspace{11mu} \xi} = {{\sum\limits_{i = 1}^{m}\; u_{ci}} + v_{ci} + u_{mi} + v_{mi} + u_{yi} + v_{yi}}} & \left( {{Eq}.\mspace{14mu} 142} \right)\end{matrix}$

subject to

$\begin{matrix}{{{\sum\limits_{j = 1}^{n + 4}\; {\left( {a_{cj} - b_{cj}} \right){F_{cj}\left( z_{i} \right)}}} + u_{ci} + v_{ci}} = C_{i}} & \left( {{Eq}.\mspace{14mu} 143} \right) \\{{{\sum\limits_{j = 1}^{n + 4}\; {\left( {a_{mj} - b_{mj}} \right){F_{mj}\left( z_{i} \right)}}} + u_{mi} + v_{mi}} = M_{i}} & \left( {{Eq}.\mspace{14mu} 144} \right) \\{{{\sum\limits_{j = 1}^{n + 4}\; {\left( {a_{yj} - b_{yj}} \right){F_{yj}\left( z_{i} \right)}}} + u_{yi} + v_{yi}} = Y_{i}} & \left( {{Eq}.\mspace{14mu} 145} \right) \\{{where}{a_{cj},b_{cj},a_{mj},b_{mj},a_{yj},b_{yj},u_{i},{v_{i} \geq {0\mspace{14mu} {and}}}}{F_{ji} = \left\lbrack {1,{z(i)},{\phi_{1}\left( z_{i} \right)},{\phi_{2}\left( z_{i} \right)},\ldots,{\phi_{n}\left( z_{i} \right)}} \right\rbrack_{n + 4}}} & \;\end{matrix}$

for the C, M, and Y components respectively.

When performing the two-stage BR algorithm, the vectors relating to C,M, and Y are generally be exchanged with the vectors of C, M, and Yrespectively. This yields an exemplary condensed initial simplex tablefor CIELab to CMY conversion as shown below.

One Dimensional Color Manifold in the CMYK Color Space

One of the more difficult and challenging tasks in color printing isreducing the cost of printing while maintaining quality. Traditionally,the printer color conversion from CMY (or CIELab) to CMYK has beenachieved (UCR) and (UCA). As discussed, UCR is the process of removingan equal amount of cyan, magenta, and yellow, and adding the same amountof black. More explicitly, the color components are defined as:

CMY={C,M, Y},CMYK={C,M, Y,K}, and Lab={L, a, b}.

The smallest value of the C, M, Y components is formed as K=min{C, M, Y}and K is used to convert a CMY value to a CMYK value as follows:

C′=C−K, M′=M−K, and Y′=Y−K.

Because of the nonlinearity of toners/inks, this UCR model is inaccurate, and introduces unacceptable color conversion errors. Colorsresulting from the UCR process are generally dull, muddy, and hueshifted. To mitigate this problem, the UCR process is followed by theUCA correction which consists of adding a relatively small amount ofcyan, magenta, and yellow back to the C′, M′, Y′ values respectively tomake colors richer and less hue shifted. However, this method isempirical and labor intensive and the conversion error is still quitelarge.

A solution to this color toner/ink reduction problem consists ofapproximating a nonlinear continuous and differentiable functionsmapping from CIELab to CMY and CMYK to CIELab. To optimize CMY to CMYKconversion, it is assumed that accurate continuous functions describingCMYK to CIELab conversion and CIELab to CMY conversion are as follows:

-   g(L,a,b)→(C,M,Y) ε    and f(C, M, Y)→(L,a,b).-   For a CIELab value within the printer gamut, Lab₀=(L₀, a₀, b₀),    g(L₀, a₀, b₀)=(C₀,M₀, Y₀).-   If C₀>0 & M₀>0 & Y₀>0, a (C′₀, M′₀, Y′₀, K′₀) value may be    determined such that

$\begin{matrix}{\underset{\underset{X}{}}{f\left( {C_{0}^{\prime},M_{0}^{\prime},Y_{0}^{\prime},K_{0}^{\prime}} \right)} \approx {\underset{\underset{Y}{}}{\left( {L_{0},a_{0},b_{0}} \right)}.}} & \;\end{matrix}$

Then, the CMYK value x_(n)=(C′_(n), M′_(n), Y′_(n), K′_(n)) may beapproximated iteratively such that, for n large enough,

∥f(x _(n))−y _(n)∥≦ε  (Eq. 146)

for some prescribed tolerance ε. Newton's method may be used if aninitial point is close to the minimizer algorithm. To determine aninitial point for Newton's method a value x₀ is determine such that

err=∥f(x ₀)−y ₀∥≦tol,   (Eq. 147)

where tol is an initial acceptable tolerance.

One of the problems in traditional GCR algorithms is that the CIELabvalue having equal amounts of cyan, magenta, and yellowtoner/ink isoften not equal to the CIELab value of that same amount of blacktoner/ink. In other words, the CIELab value of K generally correspondsto unequal amounts of C, M, Y. If the equal amount of C, M, Y, K valuesare used as the initial point to Newton's algorithm, the aboveinequality may not be satisfied. A better initial point may be obtainedby removing unequal amounts of cyan, magenta, and yellow from the CMYvalue respectively, and adding the K value. The maximum error using thisprocess may be improved from the typical UCR error value of 14.42 to10.50.

On this region, for any K, the associated CIELab value is computed using

f(0, 0, 0,K)→(L^(K), a^(K), b^(K))   (Eq. 148)

Then, the CMY value corresponding to this K is determined using

g(L^(K), a^(K), b^(K))→(C^(K), M^(K), Y^(K))   (Eq. 149)

to establish an equivalence

K

(C^(K), M^(K), Y^(K)).   (Eq. 150)

This repeated for a sequence of K values span the range of (0, 1). Forany given (C,M, Y) value, there exists a maximum K_(p) value whosecorresponding (C^(Kp),M^(Kp), Y^(Kp)) satisfies

C=C^(Kp) or M=M^(Kp) or Y=Y^(Kp)   (Eq. 151)

A sequence of K values is determined such that

K={K_(i) ≦K _(p) , i=1, 2, . . . , p}. tm (Eq. 152)

The corresponding (C^(Ki),M^(Ki), Y^(Ki))'s can also be removed from the(C, M, Y) value. Then, the CMY value (C, M, Y) is converted to asequence of CMYK values such that

CMYK′={(C ^(K′) ^(i) ,M ^(K′) ^(i) ,Y ^(K′) ^(i) ), i=1, 2, . . . , p},  (Eq. 153)

where C^(K′) ^(i) =C−C_(K) _(i) , M^(K) ^(i) =M−M_(K) _(i) , and Y^(K′)^(i) =Y−Y_(K) _(i) , which yields an optimizing GCR algorithm.

Given an initial CMY, the optimized GCR algorithm may determine anassociated (C′, M′, Y′, K′) such that the quality

∥f(C′,M′,Y′,K′)−f(C,M,Y,0)∥  (Eq. 154)

is relatively small. With an initial point close enough to theminimizer, the multidimensional Newton's method may be used to determinethe optimal value of (C′,M′,Y′,K′). For example, let

$\begin{matrix}{{f = \begin{pmatrix}f_{L} \\f_{a} \\f_{b}\end{pmatrix}},{x = \begin{pmatrix}C \\M \\Y \\K\end{pmatrix}},{y = \begin{pmatrix}L \\a \\b\end{pmatrix}},{y = \begin{pmatrix}L_{0} \\a_{0} \\b_{0}\end{pmatrix}},} & \left( {{Eq}.\mspace{14mu} 155} \right)\end{matrix}$

where y₀ is the corresponding CIELab value of the original CMY color(e.g., the desired CIELab value for the optimizer of the (C′,M′,Y′,K′))f(x) then represents the mapping from CMYK to CIELab and expresses as aRBF. For example,

$\begin{matrix}{{{f_{j}(x)} = {{A_{j}x} + a_{j\; 0} + {\sum\limits_{i = 1}^{N}\; {w_{ji}{\varphi_{j}\left( {{x - c_{i}}} \right)}}}}},} & \left( {{Eq}.\mspace{14mu} 156} \right)\end{matrix}$

where j=1, 2, 3 correspond to {f_(L), f_(a), f_(b)}. Then,f(x+h)≈f(x)+f(x)h, where f′(x) is in the form of Jacobian matrix

$\begin{matrix}{J = {{f^{\prime}(x)} = {\begin{pmatrix}\frac{\partial f_{L}}{\partial C} & \frac{\partial f_{L}}{\partial M} & \frac{\partial f_{L}}{\partial Y} & \frac{\partial f_{L}}{\partial K} \\\frac{\partial f_{a}}{\partial C} & \frac{\partial f_{a}}{\partial M} & \frac{\partial f_{a}}{\partial Y} & \frac{\partial f_{a}}{\partial K} \\\frac{\partial f_{b}}{\partial C} & \frac{\partial f_{b}}{\partial M} & \frac{\partial f_{b}}{\partial Y} & \frac{\partial f_{b}}{\partial K}\end{pmatrix}.}}} & \left( {{Eq}.\mspace{14mu} 157} \right)\end{matrix}$

For each iteration k, x_(k+1)=x_(k)−J(x_(k))¹(f(y_(k))−f(y₀)). Thisiteration stops when ∥(f(x_(n))−y₀)∥≦10⁻¹².

The dimension of the Jacobian matrix is generally 3×4. In order tochange it to an invertible Jacobian matrix. The black value of the (CM,Y,K) color remains unchanged. The Jacobian matrix in Equation 149 isthen changed to

$\begin{matrix}{J = {{f^{\prime}(x)} = \begin{pmatrix}\frac{\partial f_{L}}{\partial C} & \frac{\partial f_{L}}{\partial M} & \frac{\partial f_{L}}{\partial Y} \\\frac{\partial f_{a}}{\partial C} & \frac{\partial f_{a}}{\partial M} & \frac{\partial f_{a}}{\partial Y} \\\frac{\partial f_{b}}{\partial C} & \frac{\partial f_{b}}{\partial M} & \frac{\partial f_{b}}{\partial Y}\end{pmatrix}}} & \left( {{Eq}.\mspace{14mu} 158} \right)\end{matrix}$

Solving the Equation 149 using SVD obtains

$\begin{matrix}{x_{k + 1} = {\begin{pmatrix}C_{k + 1} \\M_{k + 1} \\Y_{k + 1} \\K\end{pmatrix}.}} & \left( {{Eq}.\mspace{14mu} 159} \right)\end{matrix}$

At each iteration, check if the C, M, Y values in the converted (C, M,Y. K) analyzer to determine whether they are in the range of [0,1],(e.g., if x_(k+1) ^(i)≧1 or x_(k+1) ^(i)≦0, where i=1, 2, 3correspond to C, M and Y component respectively). If the value of thei^(th) component is out of the range,

$\begin{matrix}{x_{k + 1}^{i} = \left\{ {\begin{matrix}{1,{{x(i)} \geq 1}} \\{0,{{x(i)} \leq 0}}\end{matrix}.} \right.} & \left( {{Eq}.\mspace{14mu} 160} \right)\end{matrix}$

The i^(th) component remains unchanged for the rest of the iterations.The column that corresponds to this component from the Jacobian matrix Jis then removed, and the constraint on the black component is relaxed byadding the partial derivatives with respect to black in the matrix J.For example, if x(2)≧1, the magenta component is removed. The Jacobianmatrix J_(k−1) becomes

$\begin{matrix}{J_{k + 1} = \begin{pmatrix}\frac{\partial f_{L}}{\partial C} & \frac{\partial f_{L}}{\partial M} & \frac{\partial f_{L}}{\partial Y} \\\frac{\partial f_{a}}{\partial C} & \frac{\partial f_{a}}{\partial M} & \frac{\partial f_{a}}{\partial Y} \\\frac{\partial f_{b}}{\partial C} & \frac{\partial f_{b}}{\partial M} & \frac{\partial f_{b}}{\partial Y}\end{pmatrix}} & \left( {{Eq}.\mspace{14mu} 161} \right) \\{and} & \; \\{x_{k + 1} = {\begin{pmatrix}C_{k + 1} \\1.0 \\Y_{k + 1} \\K_{k + 1}\end{pmatrix}.}} & \left( {{Eq}.\mspace{14mu} 162} \right)\end{matrix}$

Generally, any given CMY value of (C>0 & M>0 & Y>0) converges to a CMYKvalue with the CIELab value difference practically zero.

As described above, each black level K, (e.g., CMYK=(0, 0, 0, K) can bemapped to a CMY value of (C_(k),M_(k), Y_(k)). For any CMY value (C>0 &M>0 & y>0), there exists a set of values {K_(i), i=1, 2, . . . , p} inwhich a corresponding {C^(ki),M^(ki), Y^(ki), i=1, 2, . . . , p} may beremoved from C,M, and Y, (e.g.,C≧C^(ki) & M≧M^(ki) & Y≧Y^(ki), i=1, 2, .. , , p). Generally, for one CMY value, there exists a set of CMYKvalues {C−C^(ki),M−M^(ki), Y−Y^(ki),K_(i), i=1, 2, . . . , p} that canbe used for the initial points to the multidimensional Newton'soptimization algorithm. Thereafter, the set of values f⁻¹ (L₀, a₀, b₀)may be determined as a one-dimensional manifold based on the pre-imagetheorem.

For example, X→Y be a smooth map between manifolds. Thus, a point is aregular value of f if for all x: x ε f⁻¹(y) and the map df:T_(x)X→T_(y)Y is subjective. Here, T_(x)X and T_(y)Y are the tangentspaces of X and Y at the points x and y. Again, with f: X→Y being asmooth map, y ε Y may be a regular value of f. Accordingly, x:x ε F¹(y)is a submanifold of X with dim f¹(y)=dimX−dim Y.

For color conversion, A={set of elements of CMY K} and B={set ofelements of CIELab}. A functions=(f_(L), f_(a), f_(b)) that maps A to Bis in the form of RBFs described in Eq. 156 and f is continuous anddifferentiable. For a CIELab value y ε B, if the Jacobian matrixdescribed in Eq. 157 at y is full rank, (e.g., rank is equal to 3), theny is a regular point of f. Generally, there exists a manifold S=f⁻¹(y),S ⊂ A with the dimension of one.

These one-dimensional manifolds based on the pre-image theorem are thendetermined. For example, for any CIELab value y that is inside theprinter gamut, obtain the CMY value using g(y)→(C, M, Y). Then,determine whether min(C, M, Y)=0. If so y is not a regular point sinceits Jacobian matrix will be rank deficient. Otherwise, if min(C, M,Y)>0, the initializing algorithm is performed to obtain the sequence ofthe CMYK value CMYK′ described in Eq. 155.

The Jacobian matrix is monitored at each starting point to decide if yis a regular point. If y is a regular point, the optimized GCR algorithmis performed to obtain optimized value for each starting point, (e.g.,obtain the 1-dimensional manifold S ⊂ A, and S=f⁻¹(y)). When thereexists a level set of CMYK values S=f⁻¹(y) that correspond to the sameCIELab value, y ε B, the optimal solution for the toner/ink selectioncan be obtained for a given objective function.

The objective function can be toner/ink cost based. That is, theobjective function is to minimize the cost of printing or toner/inkcoverage. For example, min(C+M+Y+K) is the minimum amount of C+M+Y+Kobtained from the level set. The data in FIG. 5 for the original CMYvalue {0.7, 1, 1} results in the total amount of toner/ink being 270%.The minimum total amount of toners/inks in the CMYK level set is 84.2%corresponding the CMYK value of {0.066186, 0.031126, 0, 0.74461}resulting in a total toner savings of 185.8%. Because values in a CMYKlevel set correspond to the same CIELab value, the toner savings isachieved with little or no quality degradation.

Color Conversion with Ink Limitations

Another challenging task in color printing relates to the limitation ofthe amount of toner/ink coverage on a page to avoid excessivebleed-through. This task is particularly important for inkjet printing.When there exists a limit for the maximum amount of toner/ink coverageallowed per spot color printing, color conversion in highly saturatedcolor regions or dark neutral color regions are generally affected. Thecolor conversion affected by this limit includes the CIELab to CMY colorconversion. The L₁ norm approximation may provide a suitable solutionfor solving the neutral color conversion problem and the L_(∞)approximation may be suitable for color conversion in saturated colorregions. For both approaches, RBFs with linear programming may beemployed. The existence of a toner/ink limitation may be added as aconstraint in the linear programming problem formulations. For example,the color conversion with the toner/ink limitation problem may now beformulated via both L₁ and L_(∞) approximation algorithms. The L₁algorithm may be a modified BR primal algorithm with an added constraintwhile the L_(∞) algorithm may be developed based on the BP dualalgorithm to extend the three-stage algorithm to a four-stage algorithm.

To analyze ink limitation problem using the L₁ approximation, g(z) isused to represent the mapping from CIELab to CMY in the form of an RBFas:

$\begin{matrix}{{{g_{k}\left( z_{j} \right)} = {{A_{k}z_{j}} + \alpha_{k\; 0} + {\sum\limits_{i = 1}^{n}\; {\alpha_{ki}{\varphi_{k}\left( {{z_{j} - c_{i}}} \right)}}}}},} & \left( {{Eq}.\mspace{14mu} 163} \right)\end{matrix}$

where k=1, 2, 3 corresponding to {g_(C), g_(M), g_(Y)}. Then, z_(j)={L,a, b}_(j), where i=1, 2, . . . , n and j=1, 2, . . . , m. The dimensionof A_(k) is 1×3 corresponding to the dimension of the domain. Thedimension of a₀ is 1. Therefore, there are n+4 coefficients in such thatEq. 163 may be written as:

$\begin{matrix}{{{g_{k}\left( z_{j} \right)} = \underset{\underset{F_{k}{(z_{j})}}{}}{\left( {z_{j},1,{\varphi_{k}\left( {{z_{j} - c_{1}}} \right)},{\varphi_{k}\left( {{z_{j} - c_{2}}} \right)},\ldots,{\varphi_{k}\left( {{z_{j} - c_{i}}} \right)}} \right)\quad}}\mspace{565mu} \underset{\underset{w_{k}}{}}{\begin{pmatrix}A_{K}^{T} \\a_{k\; 0} \\a_{k\; 1} \\\alpha_{k\; 2} \\\ldots \\\alpha_{kn}\end{pmatrix}_{n + 4}}} & \left( {{Eq}.\mspace{14mu} 164} \right) \\{\mspace{59mu} {= {\sum\limits_{i = 1}^{n + 4}\; {w_{ki}{{F_{ki}\left( z_{i} \right)}.}}}}} & \left( {{Eq}.\mspace{14mu} 165} \right)\end{matrix}$

The L₁ error for each data point of component of C, M and Y usingnonnegative variables are defined as

$\begin{matrix}{{{g_{k}\left( z_{i} \right)} - {\sum\limits_{i = 1}^{n + 4}\; {w_{ki}{F_{ki}\left( z_{i} \right)}}}} = {u_{kj} - v_{kj}}} & \left( {{Eq}.\mspace{14mu} 166} \right)\end{matrix}$

and the coefficients for each component using nonnegative variables aredefined as w_(ki)=b_(ki)−c_(ki), where u_(kj), v_(kj), b_(ki), c_(ki)≧0.The primal linear programming problem in the L₁ sense may be illustratedby the constraint of the toner/ink limit applying to the color data witha toner/ink summation exceeding the total amount of toner or ink (e.g.,C+M+Y+K>L, where L is the total amount of toner or ink). In this regard,

$\begin{matrix}{{{Minimize}\mspace{14mu} \xi} = {{\sum\limits_{j = 1}^{m}\; u_{1\; j}} + v_{1\; j} + u_{2\; j} + v_{2\; j} + u_{3\; j} + v_{3\; j}}} & \left( {{Eq}.\mspace{14mu} 167} \right)\end{matrix}$

subject to

$\begin{matrix}{{{\sum\limits_{i = 1}^{n + 4}\; {\left( {a_{1\; j} - b_{1\; j}} \right){F_{1\; j}\left( z_{j} \right)}}} + u_{1\; j} + v_{1\; j}} = C_{j}} & \left( {{Eq}.\mspace{14mu} 168} \right) \\{{{\sum\limits_{i = 1}^{n + 4}\; {\left( {a_{2\; j} - b_{2\; j}} \right){F_{2\; j}\left( z_{j} \right)}}} + u_{2\; j} + v_{2\; j}} = M_{j}} & \left( {{Eq}.\mspace{14mu} 169} \right) \\{{{\sum\limits_{i = 1}^{n + 4}\; {\left( {a_{3\; j} - b_{3\; j}} \right){F_{3\; j}\left( z_{j} \right)}}} + u_{3\; j} + v_{3\; j}} = Y_{j}} & \left( {{Eq}.\mspace{14mu} 170} \right) \\{{{{\sum\limits_{i = 1}^{n + 4}\; {\left( {a_{1\; j} - b_{1\; j}} \right){F_{1\; j}\left( z_{jl} \right)}}} + {\sum\limits_{i = 1}^{n + 4}\; {\left( {a_{2\; j} - b_{2\; j}} \right){F_{2\; j}\left( z_{jl} \right)}}} + {\sum\limits_{i = 1}^{n + 4}\; {\left( {a_{3\; j} - b_{3\; j}} \right){F_{3\; j}\left( z_{jl} \right)}}}} \leq L},} & \left( {{Eq}.\mspace{14mu} 171} \right) \\{{{{where}\mspace{14mu} J} = {\left\{ {j_{l}{{C_{jl} + M_{jl} + Y_{jl}} > L}} \right\} \subseteq \left\{ {{i = 1},2,\ldots,m} \right\}}},u_{kj},v_{kj},b_{ki},{c_{ki} \geq 0},{k = 1},2,3,{i = 1},2,\ldots,n,{{{and}\mspace{14mu} j} = 1},2,\ldots,{m.}} & \;\end{matrix}$

With this in mind, a condensed initial simplex table for the CIELab toCMY conversion with the toner/ink limitation L is established in thetable below where p denotes the number of color data points that exceedthe toner/ink limit. The slack variables s_(i)'s are then added in Eq.171.

The slack variables si's are then added in Eq. 171.

The modified BR's algorithm may then be described in two stages. In thefirst stage, there are three sets of n iterations for the C, M, and Ycomponents respectively. In each of the n iterations, the pivotalcolumns are restricted to b_(i) and c_(i). The vector entering the basisis the one with the largest nonnegative marginal cost. The vectorleaving the basis is chosen from v_(j) and u_(j) for the one causing themaximum reduction in the objective function. The interchange of vectorsleaving and entering the basis are these performed within the same colorcomponent.

In the next stage, the nonbasic u_(j) or v_(j) are interchanged with thebasic u_(j), v_(j), or s_(j1). Neither b_(i) or c_(i) is allowed toleave the basis. The criteria for selecting vectors leaving and enteringthe basis are the same as that in the first stage. The interchange ofvectors v_(j) and u_(j) leaving and entering the basis may be performedwithin the same color component. However, w_(j1) can be interchangedwith the v_(j) and u_(j) vector for any color component. Once the vectors_(j1) leaves the basis, it is generally not allowed to enter back ontothe basis. The iteration stops when all the marginal costs arenon-positive. Thereafter, basic vectors b_(i) or c_(i) may beinterchanged with the corresponding nonbasic vectors c_(j) or b_(j).

As discussed above, the BP dual algorithm may be a relatively goodalgorithm for color conversion using L_(∞) approximation. For the colorconversion with the toner/ink limitation, a constraint may be added inthe problem formulation of the BP dual algorithm using the L_(∞)approximation.

In this regard, the L_(∞) errors may be defined for component of C, Mand Y:

$\begin{matrix}{\xi_{c} = {\max\limits_{1 \leq j \leq m}{{{\sum\limits_{i = 1}^{n + 4}\; {w_{ci}{F_{ci}\left( z_{j} \right)}}} - C_{j}}}_{\infty}}} & \left( {{Eq}.\mspace{14mu} 172} \right) \\{\xi_{m} = {\max\limits_{1 \leq j \leq m}{{{\sum\limits_{i = 1}^{n + 4}\; {w_{mi}{F_{mi}\left( z_{j} \right)}}} - M_{j}}}_{\infty}}} & \left( {{Eq}.\mspace{14mu} 173} \right) \\{\xi_{y} = {\max\limits_{1 \leq j \leq m}{{{\sum\limits_{i = 1}^{n + 4}\; {w_{yi}{F_{yi}\left( z_{j} \right)}}} - Y_{j}}}_{\infty}}} & \left( {{Eq}.\mspace{14mu} 174} \right)\end{matrix}$

Then, the problem may be formulated in the primal form by minimizingξ=ξ_(c)ξ_(m)+μ_(y), subject to

$\begin{matrix}{{\xi_{c} + {\sum\limits_{i = 1}^{n + 4}\; {\alpha_{ci}{F_{ci}\left( z_{j} \right)}}}} \geq C_{j}} & \left( {{Eq}.\mspace{14mu} 175} \right) \\{{\xi_{c} - {\sum\limits_{i = 1}^{n + 4}\; {\alpha_{ci}{F_{ci}\left( z_{j} \right)}}}} \geq {- C_{j}}} & \left( {{Eq}.\mspace{14mu} 176} \right) \\{{\xi_{m} + {\sum\limits_{i = 1}^{n + 4}\; {\alpha_{mi}{F_{mi}\left( z_{j} \right)}}}} \geq M_{j}} & \left( {{Eq}.\mspace{14mu} 177} \right) \\{{\xi_{m} - {\sum\limits_{i = 1}^{n + 4}\; {\alpha_{mi}{F_{mi}\left( z_{j} \right)}}}} \geq {- M_{j}}} & \left( {{Eq}.\mspace{14mu} 178} \right) \\{{\xi_{y} + {\sum\limits_{i = 1}^{n + 4}\; {\alpha_{yi}{F_{yi}\left( z_{j} \right)}}}} \geq Y_{j}} & \left( {{Eq}.\mspace{14mu} 179} \right) \\{{\xi_{y} - {\sum\limits_{i = 1}^{n + 4}\; {\alpha_{yi}{F_{yi}\left( z_{j} \right)}}}} \geq {- Y_{j}}} & \left( {{Eq}.\mspace{14mu} 180} \right) \\{{{- {\sum\limits_{i = 1}^{n + 4}\; \left( {{\alpha_{ci}{F_{ci}\left( z_{jl} \right)}} + {\alpha_{mi}{F_{mi}\left( z_{jl} \right)}} + {\alpha_{yi}{F_{yi}\left( z_{jl} \right)}}} \right)}} \leq {- L}},} & \left( {{Eq}.\mspace{14mu} 181} \right) \\{{{where}\mspace{14mu} \xi \; c},{\xi \; m},{{\xi \; y} \geq 0},\alpha_{1},\alpha_{2},{\ldots \mspace{11mu} \alpha_{n + 4}\mspace{14mu} {unrestricted}},{z_{j} = \left\{ {L,a,b} \right\}_{j}},{j = 1},2,\ldots,m,{J = {\left\{ {{jl}{{C_{jl} + M_{jl} + Y_{jl}} > L}} \right\} \subseteq \left\{ {{j = 1},2,\ldots,m} \right\}}},{{{and}\mspace{14mu} F_{j}} = \left\lbrack {1,z_{j},{\phi_{1}\left( z_{j} \right)},{\phi_{2}\left( z_{j} \right)},\ldots,{\phi_{n}\left( z_{j} \right)}} \right\rbrack_{n + 4}},} & \;\end{matrix}$

Transformed from the primal formulation, the dual formulation isperformed by maximizing

$\begin{matrix}{{{\sum\limits_{i = 1}^{m}\; \left\lbrack {{\left( {{\underset{\_}{\sigma}}_{ci} - {\underset{\_}{\tau}}_{ci}} \right)C_{i}} + {\left( {{\underset{\_}{\sigma}}_{mi} - {\underset{\_}{\tau}}_{mi}} \right)M_{i}} + {\left( {{\underset{\_}{\sigma}}_{yi} - {\underset{\_}{\tau}}_{yi}} \right)Y_{i}}} \right\rbrack} - {L{\sum\limits_{i_{k} = 1}^{p}\; {\underset{\_}{\rho}}_{i_{k}}}}},} & \left( {{Eq}.\mspace{14mu} 182} \right)\end{matrix}$

subject to

$\begin{matrix}{\sum\limits_{j = 1}^{m}\; \left\lbrack {{{\left( {{\underset{\_}{\sigma}}_{cj} - {\underset{\_}{\tau}}_{cj}} \right)F_{c_{ij}}^{T}} - {\sum\limits_{{jl} = 1}^{p}\; {{\underset{\_}{\rho}}_{j}F_{c_{{ij}_{1}}}^{T}}}} = 0} \right.} & \left( {{Eq}.\mspace{14mu} 183} \right) \\{\sum\limits_{j = 1}^{m}\; \left\lbrack {\left( {{\underset{\_}{\sigma}}_{cj} + {\underset{\_}{\tau}}_{cj}} \right) \leq 1} \right.} & \left( {{Eq}.\mspace{14mu} 184} \right) \\{\sum\limits_{j = 1}^{m}\; \left\lbrack {{{\left( {{\underset{\_}{\sigma}}_{mj} - {\underset{\_}{\tau}}_{mj}} \right)F_{m_{ij}}^{T}} - {\sum\limits_{{jl} = 1}^{p}\; {{\underset{\_}{\rho}}_{j}F_{m_{{ij}_{1}}}^{T}}}} = 0} \right.} & \left( {{Eq}.\mspace{14mu} 185} \right) \\{\sum\limits_{j = 1}^{m}\; \left\lbrack {\left( {{\underset{\_}{\sigma}}_{mj} + {\underset{\_}{\tau}}_{mj}} \right) \leq 1} \right.} & \left( {{Eq}.\mspace{14mu} 186} \right) \\{\sum\limits_{j = 1}^{m}\; \left\lbrack {{{{\left( {{\underset{\_}{\sigma}}_{yj} - {\underset{\_}{\tau}}_{yj}} \right)F_{y_{ij}}^{T}} - {\sum\limits_{{jl} = 1}^{p}\; {{\underset{\_}{\rho}}_{j}F_{y_{{ij}_{1}}}^{T}}}} = 0},} \right.} & \left( {{Eq}.\mspace{14mu} 187} \right) \\{\sum\limits_{j = 1}^{m}\; \left\lbrack {{\left( {{\underset{\_}{\sigma}}_{yj} + {\underset{\_}{\tau}}_{yj}} \right) \leq 1},} \right.} & \left( {{Eq}.\mspace{14mu} 188} \right) \\{and} & \; \\{{\underset{\_}{\sigma}}_{j},{\underset{\_}{\tau}}_{j},{{\underset{\_}{\rho}}_{j} \geq 0},{j = 1},2,\ldots,m,} & \left( {{Eq}.\mspace{14mu} 189} \right)\end{matrix}$

where z_(j)={L, a, b}_(j), j=1, 2, . . . , m,J={j1|C_(j1)+M_(j1)+J_(j1)>L}⊂{j=1,2, . . . , m}. F_(j)=[1, z_(j),φ₁(z_(j)), φ₂(z_(j)), . . . , φ_(n)(z_(j))]_(n+4), F_(j1)=[1, z_(j1),φ_(j1)(z_(j1)), φ_(j2)(z_(j1)), . . . , φ_(n)(z_(j1))]_(n+4).φ_(ij)=φ_(i)(z_(j)), and p denotes the number of color data points thatexceed the toner/ink limit.

A modified BP L_(∞) approximation algorithm may be used for solving thecolor conversion from CIELab to CMY with a toner/ink limitation L basedon the condensed initial simplex table below.

In the first stage of the BP algorithm, the first n simplex iterationsare performed. For example, the a_(c) vectors may be moved out of thebasis. The vector entering the basis is selected from s_(c) tocorrespond to that with the largest absolute reduced cost. The a_(m)vectors may be moved out of the basis. The vector entering the basis isselected from s_(m) to correspond to that with the largest absolutereduced cost. The a_(y) vectors may be moved out of the basis. Thevector entering the basis is selected from s_(y) to correspond to thatwith the largest absolute reduced cost.

In the second stage of the BP algorithm, ξ_(c) is forced out of thebasis. The vector entering the basis is selected from the remainingnonbasic vectors s_(c) and a_(c) is not allowed back onto the basis.ξ_(m) is forced out of the basis. The vector entering the basis isselected from the remaining nonbasic vector s_(m) and a_(m) is notallowed to move back onto the basis. ξ_(y) is forced out of the basis.The vector entering the basis is selected from the remaining nonbasicvector s_(y) and a_(y) is not allowed to move back onto the basis.

In the third stage of the BP algorithm, the pivotal column is chosencorresponding to the most negative marginal cost. The pivotal row maythen be chosen by the ratio selection rule. The interchange of nonbasicvectors s_(i) and basic vectors s_(j) may be performed within the samecolor component. The iteration continues until all marginal costscorresponding to the nonbasic vectors s_(c), s_(m), and s_(y), arenonnegative. The vector s_(cmy) is not allowed to enter the basis.

A fourth stage is added to perform simplex iteration interchangingnonbasic vectors including all the remaining vectors s_(c), s_(m),s_(y), and s_(cmy) and the basic vectors in the basis. Nonbasic vectorss_(c), s_(m), and s_(y) interchange with basic vectors of s_(c), s_(m),and s_(y) respectively. The nonbasic vectors s_(c), s_(m), and s_(y)interchange with any basic vectors of s_(c), s_(m), and s_(y) endingwhen the marginal costs are nonnegative. The function φ thus satisfiesthe Haar condition as described in above.

Gamut Mapping

The inventive embodiments of out-of-gamut conversion that include CIELabto CMYK gamut conversion may be described as follows. First, thethresholds for color difference in hue angle Th and color difference inlightness TL based on the perceptual color difference tolerance aredefined. Next, the output device gamut boundary gmt0 specified in theCIELab color space is defined. This boundary may be expanded into a fewlayers with a small increment each time (e.g., {gmt_(i), i=1, 2, . . . ,N}). In each layer, the CIELab values are {(L, a, b)_(ij), j=1, 2, . . ., m_(i)}. The corresponding color in CIELch is {(L, c, h)_(ij), j=1, 2,. . . , m_(i)} where mi is the number of the CIELab colors in the layeri. Then the CMY value is obtained using

g₀((L, a, b)_(1j))→(C,(M Y )_(1j),   (Eq. 190)

where g0 is the mapping function from CIELab to CMY based on the insidegamut colors. The CIELab data in the first expanded layer is (L, a,b)₁={(L, a, b)_(1j), j=1, 2, . . . , m₁}. The number of color datapoints in gmt1 is (C,M,Y)₁={(C,M,Y)_(1j), j=1, 2, . . . , m₁}, m₁. Ifany C, M, and Y component value in (C, M, Y)1 is smaller than 0 orgreater than 1, then the value is changed to 0 or 1 respectively. The(C, M, Y )1 may then be changed to (C, M, Y)′1. (L, a, b)′1 is obtainedusing

f₀(C′_(1j),M′_(1j),Y′_(1j),0)→(L′_(1j),a′_(1j),b′_(1j))   (Eq. 191

where f₀ is the mapping function from CMYK to CIELab based on the insidegamut colors. Then (L, a, b)′1 to (L, C, h)′1 is converted and thelightness difference ({dL_(1j), j=1, 2, . . . , m₁} and the huedifference ({dh_(1j), j=1, 2, . . . , m₁} are calculated for each datapoint between (L, C, h)₁ and (L, C, h)′₁. If dh_(1j)>T_(h) ordL1j>T_(I), (L, C, h)′₁ is converted to (L, C, h)″₁ by making thefollowing changes in L′_(1j) and h′_(1j):

$\begin{matrix}{L = \left\{ \begin{matrix}{{L_{1\; j} + T_{L}},{L_{1\; j}^{\prime} > L_{1\; j}}} \\{{L_{1\; j} - T_{L}},{L_{1\; j}^{\prime} < L_{1\; j}}}\end{matrix} \right.} & \left( {{Eq}.\mspace{14mu} 192} \right) \\{h = \left\{ {\begin{matrix}{{h_{1\; j} + T_{h}},{h_{1\; j}^{\prime} > h_{1\; j}}} \\{{h_{1\; j} - T_{h}},{h_{1\; j}^{\prime} < h_{1\; j}}}\end{matrix}.} \right.} & \left( {{Eq}.\mspace{14mu} 193} \right)\end{matrix}$

The steps may be repeated until the (C, M, Y)₁ values are found to meetthe criteria of T_(h) and T_(L). Then all out-of-gamut CIELab values (L,a, b)₁ may be mapped in gmt1 with (C, M, Y)₁. Models of g₁ and f₁ arethen generated using the data points in gmt1:

g₁((L,a,b)₁)→(CMY) and   (Eq. 194)

f₁((C,M,Y)1)→(L,a,b)   (Eq. 195)

Each of the above steps may then be repeated and the CIELab datainterpolated in {gmt_(i), i=2, 3, . . . , N} to find {(C,M,Y)_(ij), j=1,2, . . . , m_(i)}. The iteration stops when the last color layer reachesoptimized color ranges (e.g., as described in U.S. Pat. No. 6,719,392,the entire contents of which are incorporated by reference)

When the RBF models f_(i) and g_(i) are generated based on the datapoints in the gamut layer gmt_(i), the data used in generating thesemodels may contain two parts. The first part includes the data that fitsthe model g_((i−1)) and the data of f⁽¹⁻¹⁾ that is generated based onthe gamut layer i−1. The second part includes the data not fittingg_((i−1)) and f_((i−1)). This second part has C, M, or Y values eitherbelow 0 or above 1 using f_((i−1)) and is clipped to either 0 and 1 andrefined with the thresholds T_(h) and T_(L).

In doing so, the models generated for each gamut layer are continuousand smooth. To summarize, the gamut mapping direction is obtained bycreating the numerical model guided by a perceptual color differencemodel. The advantages of this technique are that more levels areobtained from this gamut mapping algorithm. Moreover, out-of-gamutcolors in a close neighborhood vary continuously. The direction of thegamut mapping is not fixed, rather it varies smoothly to provideaesthetically pleasing color mappings.

FIG. 6 is a flowchart illustrating a method 600 of optimizing colorconversion in one exemplary embodiment of the invention. In thisembodiment, the modeling module 101 receives image data that includescolor values in the process element 601. From there, the modeling module101 determines color regions in the image data in the process element602. For example, the modeling module may determine certain areas in theimage data that contain neutral color regions, medium color regions, andsaturated color regions that each requiring approximations includingvarious errors that may be introduced during color conversion. In thisregard, the modeling module 101 may generate various numerical modelsfor each of these color regions using RBFs in the process element 603.These models may be used by the converter 102 to convert the colorvalues within the regions to a desired color space, such as that of aprinter.

The converter 102 may also identify overlapping regions in the imagedata in the process element 604. To ensure that color conversions arecontinuous and smooth, the converter 102 may process the overlappingregions in the process element 605. For example, the converter 102 mayoptimize color values within those color regions by including colorvalues from adjacent regions in the numerical model used to convert aparticular color region such that the transition between the convertedcolor region and the adjacent regions is smooth and continuous.

Prior to color conversion, the validation module 103 may estimate thenumerical models' abilities to convert the color regions to the desiredcolor space in the process element 606. Thus, the validation module 103may generate color conversion estimates for the color values within eachregion and determine whether the converted colors accurately representthe color values of the image data in the process element 607. If thecolor conversion estimates do not reflect accurate representations ofthe image data color values, the validation module 103 optimizes thenumerical model(s) of the color region(s) in the process element 609.For example, the validation module 103 may optimize the numerical modelsby optimizing the color conversion estimates generated as L₁ errorapproximations for neutral color regions, L₂ error approximations formedium range color regions, and L_(∞) error approximations for saturatedcolor regions. Thereafter, the method 600 returns to the process element606 such that the validation module 103 may perform another crossvalidation of the color conversion. If the color conversion of aparticular color region is deemed accurate, the converter 102 proceedsto convert the color values at the image data within that color regionto the desired color space in the process element 608 and prepare theconverted image data for printing in the process element 610.

Although specific embodiments were described herein, the scope of theinvention is not limited to those specific embodiments. The scope of theinvention is defined by the following claims and any equivalentsthereof. For example, the context of the detailed description oftenrelates to color conversion of image data from an imaging device, suchas a camera, for the purposes of printing in a different, typicallysmaller color gamut. However, the invention is not intended to be solimited. Rather, the color conversion described herein may be used toconvert image data from virtually any color space to another.

In the context of printing, the invention is not intended to be limitedto any particular type of printer. For example, the color conversiondescribed herein may be employed in either ink or toner based printers.Moreover, the color conversion may be employed in printers ranging insize from small household printers to large commercial printers. In oneembodiment, the color conversion herein may be implemented as softwareinstructions operable with a host system to convert colors betweenvarious spaces as needed. For example, software instructions may directa processing system, such as a general purpose computer, to convertimage data from one device to the color space of another device prior topresentation to the other device.

1. A method of converting image data to a color space of an outputimaging: device, the method comprising: identifying color regions in theimage data; identifying portions of the color regions that overlap inthe image data; processing the overlapping portions of the color regionsto remove abrupt transitions between the overlapping portions of thecolor regions; and converting color values in the color regions to thecolor space of the output imaging device based on a plurality ofnumerical models.
 2. The method of claim 1, further comprisinggenerating the plurality of numerical models according to a plurality ofradial basis functions.
 3. The method of claim 2, wherein generating theplurality of numerical models according to a plurality of radial basisfunctions comprises generating an error approximation of each colorregion in the image data.
 4. The method of claim 3, wherein generatingan approximation of each color region comprises: generating L₁ errorapproximations for neutral color regions; generating L₂ errorapproximations for medium range color regions; and generating L_(∞)error approximations for saturated color regions.
 5. The method of claim1, wherein processing the overlapping portions comprises optimizing theoverlapping portions by converting color values of one color region tothe color space using color values of an adjacent region.
 6. The methodof claim 1, further comprising validating color conversion of thenumerical models prior to converting the color values to select colorconversion algorithms for the color regions.
 7. The method of claim 6,wherein validating color conversion of the numerical models comprisesoptimizing color conversion estimates generated as L₁ errorapproximations for neutral color regions, L₂ error approximations formedium range color regions, and L_(∞) error approximations for saturatedcolor regions.
 8. The method of claim 1, further comprising generatingthe plurality of numerical models according to human perception and asource of printer errors.
 9. A color conversion system, comprising: amodeling module operable to identify color regions in image data anddetermine portions of the color regions that overlap, wherein themodeling module is further operable to generate a plurality of numericalmodels used in converting the color regions to a color space of anoutput imaging device; and a converter operable to convert color valuesin the color regions to the color space of the output imaging devicebased on the numerical models, wherein color conversion of the colorvalues of one color region is performed according to a color conversionalgorithm for said one color region that includes color values from anoverlapping color region.
 10. The color conversion system of claim 9,wherein the numerical models include radial basis functions.
 11. Thecolor conversion system of claim 9, further comprising a validationmodule operable to estimate color conversion of the numerical modelsprior to color conversion by the converter.
 12. The color conversionsystem of claim 11, wherein the validation module is further operable tooptimize color conversion estimates generated as L₁ error approximationsfor neutral color regions, L₂ error approximations for medium rangecolor regions, and L_(∞) error approximations for saturated colorregions.
 13. The color conversion system of claim 9, wherein themodeling module is further operable to generate L₁ error approximationsfor neutral color regions, L₂ error approximations for medium rangecolor regions, and L_(∞) error approximations for saturated colorregions.
 14. The color conversion system of claim 9, wherein themodeling module is further operable to generate the plurality ofnumerical models according to human perception and a source of printererrors.
 15. A printing system, comprising: a printer; and a printercontroller operable to process image data having a plurality of colorvalues and convert the color values of the image data to a color spaceof the printer, wherein the printer controller is further operable toidentify color regions in the image data, determine portions of thecolor regions that overlap, generate a plurality of numerical modelsused in converting the color regions to the color space of the printer,and convert the color values of the color regions to the color spacebased on the numerical models for printing to the printer, wherein colorconversion of the color values of one color region is performedaccording to a color conversion algorithm for said one color region thatincludes color values from an overlapping color region.
 16. The printingsystem of claim 15, wherein the numerical models include radial basisfunctions.
 17. The printing system of claim 15, further comprising avalidation module operable to estimate color conversion of the numericalmodels prior to color conversion by the converter.
 18. The printingsystem of claim 17, wherein the validation module is further operable tooptimize color conversion estimates generated as L₁ error approximationsfor neutral color regions, L₂ error approximations for medium rangecolor regions, and L_(∞) error approximations for saturated colorregions.
 19. The printing system of claim 15, wherein the modelingmodule is further operable to generate L₁ error approximations forneutral color regions, L₂ error approximations for medium range colorregions, and L_(∞) error approximations for saturated color regions. 20.The printing system of claim 17, wherein the modeling module is furtheroperable to generate the plurality of numerical models according tohuman perception and a source of printer errors.